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PREFACE 
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The reader is assumed to be a TSS/360 
user with a particular concern for perfor- 
mance information associated with individu- 
al subroutines. The numerical analyst is 
provided with a description of the 
algorithms, and a specification of accuracy 
and range, where these are considered to be 
significant. 
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INTRODUCTION 



The PL/I Library computational subrou- 
tines provide support for the operators and 
built-in functions of the PL/I language in 
four major categories: 

1. Bit and Character Strings 

2. Arithmetic 

3. Mathematical 

^. Arrays 

This publication gives detailed informa- 
tion in each of the four sections mentioned 
above with respect to accuracy ^ choice of 
algorithm i and range of values handled 
Cwhere appropriate) • 

A number of exceptional conditions may 
arise in the execution of the library sub- 
routines- Many of these are not directly 
related to PL/I ON conditions. The method 
of treatment in these cases is to write a 
diagnostic message and raise the ERROR con- 
aition. This allows the user the oppor- 
tunity to investigate the error by use of 
the ONCODE built-in function in his ON 



ERROR unit and to program the action he 
wants taken. 

Module Names 

The module name for each of these sub- 
routines is IHEWxxx, where xxx is usually a 
mnemonic group indicating the module 

function. 

The seventh character usually defines 
the base, scale, mode and precision of the 
arguments for a given module. In the 
arithmetic, mathematical and array subrou- 
tines, this suffix is usually one of the 
characters shown in Figure 1| the only 
exceptions to this are the array indexing 
subroutines, where the suffixes are mnemon- 
ic only, and the ALLCx) , ANY(x) subrou- 
tines, where the suffix is 1 or 2. 

In the string subroutines, the seventh 
character in each module name has only a 
ffneroonic significance. In some cases the 
seventh character may be one of those given 
in Figure 1. This is purely coincidental? 
the meanings in Figure 1 do not apply to 
the string subroutines. 



Seventh 






Character 


L 


Argument Attributes 


B 


r 
Real 


fixed-point binary 


D 


Real 


fixed-point decimal 


F 


Real 


fixed-point binary or decimal 


G 


Real 


or complex short floating-point 


H 


Real 


or complex long floating-point 


L 


Real 


long floating-point 


S 


Real 


short floating-point 


U 


Complex fixed-point binary 


V 


Complex fixed- point decimal 


X 


Complex fixed- point binary or decimal 


w 


CompJ 


.ex short floating-point 


% 


Coropl 


.ex long floating-point 



+~ 



Argument (or element of 


Maximum 


argument) Passed in 


i Precision 




1- 


Fullword 


31 


Up to 8 bytes 


15 


Binary: fullword 


Binary: 31 


Decimal: up to 8 bytes 


1 Decimal: 15 


Real: 1 fullword 


Binary: 21 


Complex: 2 full words 


1 Decimal: 6 


Real: 1 doubleword 


Binary: 53 


Complex: 2 doublewords 


[Decimal: 16 


Doubleword 


Binary: 53 




1 Decimal: 16 


Fullword 


Binary: 21 




[Decimal: 6 


2 fullwords 


, 31 


Op to 16 bytes 


i 15 


Binary: 2 fullwords 


[Binary: 31 


Decimal: Up to 16 bytes 


Decimal: 15 


2 fullwords 


Binary: 21 




Decimal: 6 


2 doublewords 


Binary: 53 


1 


Decimal: 16 



i X- 



.. . — — .-. . . .^ X. 



. ^ - A 



figure 1. Interpretation of Seventh Character in Module Names 
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CHAFFER 1; STRING OPERATIONS AND MIILT-IN FOWCTIOHS 



The liteary string package contains 
modules for handling bit and character 
string operations. Generally, a string 
function or operator is supported by only 
one module f but in the interests of effi- 
ciency some of the bit string operators are 
provided mth additional modules to deal 
with byte- aligned input data. 

A complete list of the modules provided 
in the Library string package is given in 
Figure 2. 



BIT STRING OPERATIONS 

The 'And* Operator CS) CBit Strings! 

Module Name: IHEWBSA 



Entry Point: 
Function: 



IHEBSAO 



To implement the "and* operator between 
two byte- aligned bit strings, placing the 
result in a byte-aligned target field. 

Method: 

The current length of the target string 
is set equal either to the maximum length 
of the operands, or to the maximum length 
of the target field (when truncation is 
necessary to avoid exceeding the length 
of this field). The strings are 'and'ed 
together for a length equal to the mini- 
mum of the lengths of the operands, and 
the result is extended with zeros, if 
necessary, up to the current length cal- 
culated for the target field. 



Ihe 'Or' Operator CI Halt Strings) 
Module Name: IHEUBSO 

Entry Point: IHEBSOO 

Function: 

To impleoient the 'or' operation between 
two byte-aligned bit strings, placing the 
result in a byte-aligned target field- 
Method: 

The current length of the target string 
is set equal to either the maximum length 
of the operands or to the maxiraum length 
of the target field (when truncation is 
necessary to aiioid exceeding the length 
of this field). The strings are 'or'ed 
together for a length equal to the mini- 
mum of the lengths of the operands and 
the remainder of the longer string is 
moved into the target field up to the 
current lengthi the remainder of the tar- 
get field is left unchanged. 

The 'Not' Operator 1%} (Bit Strings) 

Module Name: IHEWBSN 

Entry Point: IHEBSNO 

Function: 

To implement the 'not' operator for a 
byte-aligned bit string, placing the 
result in a byte-aligned target field. 



1 PL/I 
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1 Operation 
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Function 
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1 


String | 
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1 


General 


1 Byte -aligned! 


1 'And' (S) 
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1 


Ose BODL 
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IBJ&mSA 


1 


1 
1 


1 'Or' C|) 






« 




Use BOOL 




IHEHBSO 




1 


1 'Nof C^) 






- 




Ose BOOL 




IHEKTBSN 




1 


1 Concatenate 


MP 


1 


REPEAT 




IHEWBSR 




- 




IHEiiCSK 1 


1 Compare 






- 




IHEWBSD 




XHEWBSC 




IHIHCSC 1 


1 Assign 






- 




IHEMBSK 




IHENBSM 




IHEUCSM 1 


1 Fill 
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IHEWBSM 




- 




IHEWCSM 1 
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IHEUCSM 1 
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SOBSTR 




IHEWBSS 




- 
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INDEX 
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IHEWCSI 1 


1 -,.. ,™. - 
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Figure 2. Bit and Character String Operations and Functions 



Method: 

The current length of the target string 
is set equal to either the current length 
of the operand or to the maximum length 
of the target field Cwhen truncation is 
necessary to aiioid exceeding the length 
of this field). The target field is set 
to a string of I's for a length equal to 
its calculated current length and the 
result is obtained by an 'exclusive or* 
with the operand. The remainder of the 
target field beyond the calculated cur- 
rent length is left tmchanged. 



Concatenate/REPEM* /General Assign CBit 
Strings) 



Module Name: 



Entry Points: 



IHEieSK 



Operation 
Concatenate C | | ) 
REPEAT (Bit string, n) 
General assign 

Function: 



Entry 

Point 

IHEBSKK 

IHEBSKR 

IHEBSKA 



IHEBSKK: to concatenate two bit strings 
into a target field. 

IHEBSKR: to concatenate n ♦ 1 instances 
of the single source string into a tar- 
get field. If n < 0, the result is the 
string itself- 

IHEBSKA: to assign a bit string to a 
target field without zero filling. 

Method: 

The current length of the target field is 
made equal to the smaller of two values: 

- the sum of the current lengths of the 
source strings 

- the maximiwi length of the target field 

All entry points use a subroutine that 
cA> tains data from a source, aligns it 
correctly and moves it to the target 
field: 



IHEIKKA: Oses the subroutine once to 
mDve the source string to the target 
field. 

For all entry points, the remainder of 

the target field beyond the calculated cur- 
rent length is left unaltered. 



Comparison CBit Strinqrs, Byte-aliiynedl 
Module Name: IHEMBSC 



Entry Pointi 
Function: 



IHEBSCO 



To compare two byte- aligned bit strings 
ar^ to return a condition code as bits 2 

and 3 of a fullword target field as 

follows: 

00 if strings are equal 

01 if first string compares low at the 

first inequality 
10 if first string compares high at the 
first inequality 

The shorter string is treated as though 
extended with zeros to the length of the 
longer. 

The first byte of the target field is 
also used to preserve the program mask in 
the PSW for the calling routine. This 
byte contains: 



Bits 
to 1 
2 to 3 
*l to 7 



Contents 
Instruction length code 01 
Condition code as above 
Program mask (calling routine) 



Method: 



The two strings are compared up to the 
current length of the shorter striiig. 
The remainder of the longer string is 
cougar ed with zeros. 

General Comparison CBit Strings) 

Module Name: IHEWBSD 

Entry Point: IHEBSDO 

Function: 



IHEBSKK: uses this subroutine tmce to 
move the source strings to the target 
field. 

IHEBSKR: uses the subroutine to conca- 
tenate the contents of the target field 
with itself (whenever possible) as well 
as concatenating the contents of this 
field with the source string. Direct 
concatenation of the scwrce string 
n ♦ 1 times is not used. 



To compare two bit strings and return a 
condition code as bits 2 and 3 of a full- 
word target field as follows: 

00 if strings are equal 

01 if first string compares low at the 

first inequality 

10 if first string compares high at the 
first ine^ality 
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The shorter string is treated as though 
extended with zeros to the length of the 
longer . 

The first byte of the target field is 
also used to preserve the program mask in 
the PSW for the calling routine. This 
byte contains: 

Bits Contents 

to 1 Instruction length code 01 

2 to 3 Condition code as above 

U to 7 Prograin mask (calling routine) 



Method: 

The two strings are compared up to the 
current length of the shorter string. 
The remainder of the longer string is 
compared with zeros. 

Assign/Fill (Bit Strings) 

Module Name: IHEWBSM 



Entry Points: 



Operation 
Fixed-length assign 
Variable- length assign 
Zero fill only 

Function: 



Entry 

Point 
IHEBSMF 
IHEBSMV 
IHEBSMZ 



IHEBSMF: to assign a byte-aligned string 
to a byte-aligned fixed-length target, 
filling out with zero bits if 
necessary. 

IHEBSMV: to assign a byte-aligned string 
to a byte- aligned variable -length 
target. 

IHEBSMZ: to fill out the target area 
from its current length to its maximum 
length with zero bits. 

Method: 

IHEBSMF: the ittinirauiB of the source cur- 
rent length and the target roaximuim 
length is calculated and the source 
string is moved to the target for a 
length equal to this length. Zero 
filling of the target is performed if 
necessary. The current length of the 
target is set equal to the maximmn 
length. 

IHEBSMV: the source string is moved to 
the target field as above, but without 
zero filling. The current length of 
the target is set appropriately. 

IHEBSMZ: zeros are propagated in the 
target from the current length to the 



maximum length. The current length of 
the target is set equal to the roaucimuiii 
length. 

Otlier Information: 

This routine supplies assignment of byte- 
aligned bit strings of both fixed and 
variable lengths. Non-aligned strings 
may Jae assigned by using the general 
assign module (entry point IHEBSKA) . Any 
filling required for fixed length strings 
can then be obtained using the IHEBSMZ 
entry described above. 



BIT STRING FUNCTIONS 
SUBSTR (Bit Strings) 
Module Name: IHEWBSS 



Entry Points: 



Operation 
SUBSTR (Bit-string , i) 
SUBSTR (Bit-string, i, j) 

Puncti on : 



Entry 

Point 

IHEBSS2 

IHEBSS3 



To produce a string dope vector describ- 
ing the Sl/BSTR pseudo- variable and func- 
tion of a bit-string. 

Method: 

Arithmetic is performed according to the 
function definition, using the current 
length of the argument string. The 
result describes a fixed-length string. 

Error and Exceptional Conditions: 

STRINGRANGE 

INDEX (Bit Strings) 

Module Name: IHEWBSI 

Entry Point: IHEBSIO 

Function: 

To compare two bit strings to see if the 
second is identical to a substring of the 
first, and, if it is, to produce a binary 
integer (the index) %diich indicates the 
first bit position in the first string at 
which such a substring begins. If no 
such index is found, or if either string 
is null, the function value is zero. 

Method: 

The index is found by shifting and com- 
paring portions of the two strings in 
registers. 



BOOL CBoolean Foncrtion) CBit Strings) 

Module Name: IHEWBSF 

Entry Point: IHEBSFO 

Function: 

To take two source strings and perform 
one of the sixteen possible logical 
operations between corresponding bits. 
The particular operation performed is 
defined by inserting the bit 'pattern - 
nj^n^nan - yielded by the third argument 
into the table below: 



I First field | | | 1 | 1 | 

j.^, ^ ^^^^.„^^._^^.„^^. — -i 

I Second field | | 1 | | 1 | 

I ^ ^+^^^^^^, 4„^ 

I Target field ( n^f na| nal n,* | 

i ^, X. -JL- .^ JL _J. J 

Method: 

The current length of the target string 
is set equal to either the maximuii of the 
current lengths of the source strings or 
to the maxiwum length of the target field 
(when truncation is necessary to avoid 
ecceeding the length of this field). The 
necessary operation is performed on the 
strings and the result stored in the tar- 
get field. If one string is shorter than 
the other, it is regarded as Deing 
extended on the right with 2erc» up to 
the length of the longer. The field 
between the calculated current length and 
the maximum length of the target is left 
unchanged. 



CHARACTER STRING OPERATIONS 

Concatenate/REPEAT CCharacrter Strings) 
Module Name: IHEWCSK 



Entry 

Point 

IHECSKK 

IHECSKR 



Entry Points: 



Operation 
CoiK^atenate C||) 

REPEAT 

(Character string, n) 

Function: 



IHECSKK: to concatenate two character 
strings into a target field. 

IHECSKR: to concatenate n ■*- 1 instances 
of the single source string into a tar- 
get field. If n < 0, the result is the 
string itself. 



Method: 

The current length of the target field is 
made equal to the ana Her of two values: 

- the sum of the current lei^ths of the 

source fields. 

- the maximum length of the target 

field. 

Both entry points use a subroutine that 
moves characters from a source to the 

target: 

IHECSKK: Dses the subroutine to perform 
the required number of source moves. 

IHECSKR: Itees the subroutine to conca- 
tenate the source string with the tar- 
get field and also to concatenate the 

target field with itself (whenever 

possible) . 

For both entry points, characters beyond 
the range of the target current length 
renBin unaltered. 



Compare (Character Strings) 
Module Name: IHEWCSC 



Entry Point: 
Function: 



IHECSCO 



To compare two character strings and to 
return a condition code as bits 2 and 3 
of a fullword target field as follows: 

00 if strings are equal 

01 if first string compares low at the 

first inequality 
10 if the first string compares high at 
the first inequality 

The shorter string is treated as though 
extended with blanks to the length of the 
longer one. 

The first byte of the target field is 
also used to pr^erve the program mask in 
the PSW for the calling routine. Ttiis 
byte contains : 



Bits 
to 
2 to 
4 to 



Contents 
1 Instruction length code 01 
3 Condition code as above 
7 Program mask (calling routine) 



Method: 



Ifie two strifes are compared in storage. 

If the strings are of different lengths 
and are identical up to the length of the 
shorter, the remainier of the longer is 

compared with blanks . 
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Assiqn/Fiil/HIGH/LQW CCharacter Strings) 
Module Name: IHEWCSM 
Entry Points: 





Entry 


Operation 


Point 


Fixed-length assign 


IHECSMF 


Variable -length assign 


IHECSMV 


Blank fill only 


IflECSMB 


HIGH 


IHECSMH 


LOW 


IHECSJML 



Function: 

IHECSMF: to assign a character string to 
a fixed- length target^ filling out with 
blanks if necessary. 

IHECSMV: to assign a character string to 
a variable-l«igth target. 

IHECSMB: to fill out the target field 
from its current length to its maximum 
length with blanks. 

IHH:smB2 to fill a target field with the 
highest character in the collating 
segilience» tip to its current length* 

IHECSML: to fill the target field with 
the lowest character in the collating 
sequence, up to its current length. 

Method: 

IHECSMF: The minimum of the scwirce cur- 
rent length and the target maximum 
length is calculated and the scarce 
string is moved to the target for a 
length equal to this length. Filling 
of the target with blanks v^ to the 
target maximum length is performed if 
necessary. The current length of the 
target is set equal to its maximum 
length. 

IHECSMV: moves the string as above, but 
without blank filling. The current 
length of the target is set 
a ppr opri ate ly . 

IHECSMB: propagates blanks and sets the 
current length of the target equal to 
its maximum length - 

IHECSMH, IHECSML: uses part of the blank 
fill routine to propagate the highest 

or lowest character in the collating 
sequence up to the current length of 
the target. 



CHARACTER STRING FDNCTIONS 

SUBSTR C Character Strings) 
Module Name: IHEWCSS 



Operation 
SUBSTRCCharacter-string,il 
SOBSTR (Character-string, i, j> 

Functi on: 



Entry 

Point 

IHECSS2 

IHECSS3 



To produce a string dope vector describ- 
ing the SUBSTR pseudo- variable and func- 
tion of a character string. 

Method: 

Arithmetic is performed according to the 
function definition, using the current 
length of the argument string. The 
result describes a fixed-length string. 

Error and Exceptional Conditions: 

SIJIINGRAIKSE 



INDEX (Character Strings) 
Module Name: IHEWCSI 



Entry Point: 
Function: 



IHECSIO 



To compare two character strings to see 
if the second is identical to a substring 
of the first, and, if it is, to produce a 
binairy integer (the index) which indic- 
ates the first character position in the 
first string at which such a substring 
begins. If no siK:h index is found, or if 
either string is null, the function value 
is zero. 

Method: 

The first string is scanned from left to 
right for a character equal to the first 
character in the second string- If a 
match is found, the whole of the second 
string is compared with a substring of 
the first string beginning at the match- 
ing character. If they are equal, an 
index is produced. The scanning con- 
tinues until either an index is produced 
or the end of the fiist string is 
reached. 



CHAPTER 2: ARITHIHa^IC OPERM'IOWS AHD BUILT- HI FUHCTIOWS 



Library arithmetic modules support all 
those arithmetic generic fmctions and 
operators for i^ich the compilers neither 
produce in-line code nor Cas for the func- 
tions FIXED, FLOAT, BINARY and DECB4AL> use 
the conversion package. The names of the 
library modules %ihich support the arith- 
metic operations are given in Figure 3; the 
names of those which support the arithmetic 
functions are given in Figure <♦. 

Statistics for accuracy of floating- 
point modules are given where considered 
meaningful and helpful; an es^lanation of 
their use is given in the chapter on mathe- 
matical routines. Precise results are 
obtained from all fixed- point modules 
except complex division and complex ABS, 
idftiere small truncation errors inevitably 
occur, and the ADD function Cf ixed deci- 
mal) , in %*iich the effect of truncation 
errors depends on the relative values of 
the scale factors of the arguments. 

Any restrictions on the adtaiissibility of 
arguments are noted under the headings 
•Range* and 'Error and Exceptional 
Conditions* . 

Range z This states any ranges of arguments 
for which a modiile is valid. Arguments 
outside the ranges given are asstnned to 
have been excluded before the module is 
called. 

Error ai^l Exceptional Conditions ; These 
cover conditions which may result from the 



use of a routine; they are listed in four 

categories : 



p — Programmed conditions in the module 
concerned. Programmed tests are 
made where this is not too costly 
and, if an invalid argument is 
found, a branch is taken to the 
entry point IHEERRC of the execution 
error package CEXEPK !niis results 
in the printing of an appropriate 
message and in the ERROR condition 
beii^ raised. 

I — Interruption conditions in the 

iM>dule concerned. For those rou- 
tines where SIZE and FIXEDOVERFLOW 
are detected by programmed tests or 
where hardware interruptions may 
occur, the OVERFLOW, UNDERFLOW, 
FIXEDOVERFI/)W, SIZE and ZERODIVIDE 
conditions pass to the ON handler 
flHEERR) and are treated in the 
normal way. The machine is assumed 
to be enabled for all interruptions 
except significance, which is masked 
off. 

O ~ Programmed conditions in modules 
called by the module concerned. 
These occur when invalid arguments 
are detected in the iwDdule called. 

H — As I, but the interruption condi- 
tions occur in the modules called by 
the module concerned - 



ARITHMETIC OPERATIONS 



-^ . ^ — , J- . 

I Binary | Decimal | Short 
I fixed I fixed [ float 



^ 

Long I 
float I 

4 

I 

-—4 

IHEWXIL I 
IHEWXXLI 

I 



Operation 



Real Operations 



-^ ^ . — ^ — ^- 

I IHEWXIBJ IHEWXIDJ 

I - I ~ I 
I - I IHEWAPDl 



Integer exponentiation: x**n 
General exponentiation: x**y 
Shift-and-assign, Shift-and-load 



IHEWXIS I 

IHEWXXSI 

I 



— ^ 

I 

— 4 



h 



Complex Operations 



Multiplication/division: Zi.*2af Zj./2a I IHEWieZU| IHEWMZV| 

Multiplication: z^^z^ I "" I " I 

Division: Zj./Za I ~ I I 

Integer exponentiation: z**n | IHEWXIUJ IMEWXIV| 

General exponentiation: Zj_**Za I "" | - | 

L . , , , ^ ... . «.i , J. .±, 

Figure 3. Arithmetic Operations 



I 

IHEWMZWj 

IHEWDZWI 
IHEWXIWI 
IHEWXXWI 

. . .- X- 



IHEWMZZ I 
IHEWDZZI 
IHEWXIZ I 
IHEWXXZ I 

. . J 
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lyilTHMETIC FUNCTIONS 



I 

|. — . ^. ^ __ — ^ __ ^ _ 

t Function I Binary | Decimal | Short | Long 
I I fixed I fixed | float [float 

|. X_ i ^^ JL__ X ^ 

I Real Arguinents 

I- — 



^ 



^ ^ , — , — ^ — , ^ , — 

MAX, MIN I IHEWMXB | IKEWKXD | IHEWFJCS j IHEWMXL 
ADD I - I IHEWADD | - | - 

X X _»- ^X- ^_J 



h 



Complex Arg omenta 

^ ^_^ ^ — . .^ .^, 

I ADD I - JIHEWADVI -| - 
|MULTIPL¥|IHEWMPuilHEWFiPV| -| - 
I DIVIDE |IHEWDVU|IHEWDW| -| - 
I ABS I IHEWABU | IHEWABV | IHEWABW | IHEWAB2 

L ^_X ^ _X- -X^-^ X . .— 

Figure H^ Arithmetic Functions 



REAL OPERATIONS 

Positive Integer Exponentiation (fixed 
binary) 



Module Names IHEWXIB 
Entry Point: IHEXIBO 



Method: 

The result is set initially to the value 
of the argument. The final result is 
then obtained by repeated squaring of 
this value or squaring and multiplying by 
the argument. 

Range: 

The precision rules of PL/I impose the 
restriction that if x has a precision 
^P#q^ f this module will be called only if 
n*(p ♦ 1) - 1 < 15, This implies that 
n < 16/ fp + 1) < 8 for all such cases 
and, in fact, this module will operate 
only for the range < n < 8, 



Integer Exponentiation C floating-point) 
Module Names and Entry Points: 



Argument 
Short float 
Long float 

Function: 



Module 

Name 
IHEWXIS 
IHEWXIL 



Entry 

Point 

IHEXISO 

IHEXILO 



To calculate x**n, where n is an integer 
between -2**31 and 2 ♦♦31 - 1 inclusive. 



Functions 

To calculate x**n, where n is a positive 
integer. 

Method: 

The result is set initially to the value 
of the argument. The final result is 
then obtained by repeated squaring of 
this value or squaring and multiplying by 
the argument. 

Range: 

< n < 2^^31 

The precision rules of PL/I impose a 
further restriction in that if x has a 
precision Cp,q), this module will be 
called only if n*Cp + 1) - 1 < 31. This 
implies that n < 32/ Cp ♦ 1) < 16 for all 
such cases. 

Positive Integer Exponentiation (fixed 
decimal) 



Method : 

If the exponent is zero and the argument 
nonzero, the result 1 is returned immedi- 
ately. Otherwise the result is set ini- 
tially to the value of the argument and 
the exponent is made positive. The argu- 
ment is raised to this positive power by 
repeated squaring of the contents of the 
result field or squaring and multiplying 
by the argument. Then, if the exponent 
is negative, the reciprocal of the result 
is taken, otherwise it is left unchanged. 

Accuracy: 

The values given here are for the rela- 
tive error divided by the exponent for 
exponents between 2 and 1023; the argu- 
ments are uniformly distributed over the 
full range for each exponent for which 
neither OVERFLOW nor UNDERFLOW occurs. 
There are 2^^(10 - k) arguments for each 
exponent in the range 2*^k < exponent < 
2^^(k ♦!)-!, where k has integral 
values from 1 to 9 inclusive. 



Module Name: IHEWXID 

Entry Point: IHEXIDO 

Function: 

To calculate x^*n, where n is a positive 
integer. 



IHEWXIS 



j R.M.S. relative 
I error/exponent 
I ♦10«*6 



0.00871 



-+- 



Maximum relative | 
error/exponent | 
♦10^^6 I 



0.692 



I 



IHEWXIL 

I R.l^.S. relative | Maxiromn relatxire | 
I error/exponent j error/exponent | 
I ♦10**15 I ♦lO^^lS I 

|. ™„. — ^_^^^^^^^^^^„,„„„^ ^^4 

I 0.0995 I 1.73 I 

L . .- .^x^ .».— - — —^ . — J 



O : a. X < with y ^ 
LOG routine 



error caused in 



b. y*LOG(x) > 17^.673: error caused 
in EXP routine 

Shift-- and- assign, Shift-and-load (fixed 
deciroal) 



Error and Exceptional Conditions: 

P : X = with n < 

I : OVERFLOW, UNDERFLOW 

Since x**C-in), where m is a positive 
integer, is evaluated as l/(x**in), 
the OVERFLOW condition may occur when 
m is large, and the UNDERFLOW condi- 
tion when x is very small. 

Other Information: 

IHEWXIS: For large exponents, for 

example, those greater than 1023, it is 
generally faster and more accurate to 
use the module IHEWXXS rather than 
IHEWXIS, passing the exponent as a 
f loating-point argument. However, it 
should be noted that IHEWXXS will not 
accept a negative first argument, and 
thus it is necessary to pass the abso- 
lute value of this argument, and also, 
in cases where the exponent is odd, to 
test the sign of the argument in order 
to be able to attach the correct sign 
to the numerical result returned. 



General Floating-Point Exponentiation 



Module Names and Entry Points: 



Argument 
Short float 
Long float 

Function: 



Module 

Name 
IHEWXXS 
IHEWXXL 



Entry 

Point 

IHEXXSO 

IHEXXLO 



To calculate x**y, where x and *y are 
floating-point numbers. 

Method: 

When X = 0, the result x**y = is given 
if y > 0, and an error message if y < 0. 
When X * and y = 0, the result x**y = 1 
is given. Otherwise x**y is computed as 
EXP Cy*LOGCx) ) , using the appropriate 
mathematical function routines. 



Error and Exceptional Conditions: 

P : X = with y < 



Module Name: IHEWAPD 

Entry Points: 



Operation 
Shift and assign 
Shift and load 



Entry 

Point 

IHEAPDA 

IHEAPDB 



Function : 

IHEI^PDA: To convert a real fixed deciroal 
number with precision Cpi,gj.> to preci- 
sion Cpatga) f where pj_ < 31 and 
Pa < 15. 

IHEAPDB: To convert a real fixed deciroal 
number with precision Cpj.,gi^ to preci- 
sion (31,q2)# where p^ < 31. 

Method : 

The argument scale factor is subtracted 
from the target scale factor. The argu- 
ment is converted to precision 31 in a 
field with a shift equal to the magnitude 
of the difference between the scale fac- 
tors; the shift is to the left if the 
difference is positive and to the right 
if negative. 

If entry point IHEAPDB is used, the field 
is moved unchanged to the target. If 
entry point IHEAPDA is used, the result 
is checked for FIXEDOVERFLOW and then 
assigned to the target with the specified 
precision. The assignment may cause the 
SIZE condition to be raised. 



Error and Exceptional Conditions: 
I : FIXEDOVERFLOW or SIZE 

COMPLEX OPERATIONS 

Multiplication/Division (fixed binary) 
Module Name: IHEWMZU 
Entry Points: 



Mathematical 
Operation 

ZJ./Z2 



Entry 

Point 

IHEMZUM 

IHEMZUD 
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Fiincrtioii: 



Other Inforoiation: 



To calculate 2j.*Za or Zj^/Zat where Zi. and 
2 2 are fixed- point binary complex 
numbers . 

Method: 

Let Zi = a ♦ bl and Za = c -•- dl. Thent 
for multiplication, an incorporated sub- 
routine is used to ccanpite a*c - b*d and 
b*c ♦ a*d; these are tested for FIXED- 
OVERFLOW and then stored as the real and 
imaginary parts of the result. 

For division, the subroutine is used to 
compute a*c ♦ b*d and b*c - a*d. The 
expression c**2 ♦ d**2 is computed and 
the real and imaginary parts of the 
result are then obtained by division. 

Hie subroutine computes the expressions 

u*x ♦ v*y and v*x - u*y. 

Error and Exceptional Conditions: 

I : FIXEDOVERFLOW in either routine, 

ZERODIVIDE in the division routine* 



Multiplication/Division Cfixed decimal) 
Module Name: IHEWMEV 



Entry Points: 

Ma thema ti ca 1 
Operation 

2j./Za 



Function: 



Entry 

Point 

IHEMZVM 

IHEMZVD 



To calculate Zi.*2a or Zj./Za where Zj. 
Za are fixed-point decimal complex 
nuntoers. 



and 



Method: 



Let Zj 



a ♦ bl and Za = c -i- dl. 



The 



products a*c, b*c, a*d and b*d are com- 
puted- Then the required result is 
obtained as follows: 

Multiplication: 

Real part a*c - b*d 
Imaginary part b*c -•• a*d 

Division: 

Real part (a*c 4^ b*d)/ Cc*c ♦ d*d> 

Imaginary part Cb*c - a*d)/Cc*c ♦ d*d) 

Error and Exceptional Conditions: 

I : FIXEDOVERFLOW in either routine, 

ZERODIVIDE in the division routine. 



Where the operands differ in precision, 
it is faster to present the longer 
operand as the second argument rather 
than the first. 



Multi pli cati on ( float i nq- point ) 
Module Names and Entry Points: 



Arcpiroent 



Module 

Ha me 



Entry 
Point 



Short float IHEWMZW IHEMZWO 
Long float IHEWMZZ IHEMZZO 

Function: 

To ccHopute Zi*Za in floating-point, when 

2^ = a ♦ bl and Za = c ♦ dl. 

Method: 

The real and imaginary parts of the 
result are ccwtiputed as a*c - b*d and b* 
c ♦ a*d, respectively. 

Error and Exceptional Conditions: 

I : Exponent OVERFLOW and UNDERFLOW 

Division (floating-point) 
Module Names and Entry Points: 



Argument 



Module 
Name 



Entry 
Point 



Short float IHEWDZW IHEDZWO 
Long float IHEWDZZ IHEDZZO 

Function: 

To compute z^/Za in floating-point, when 

Zj. = a ♦ bl and Za = c ♦ dl. 

Method: 

1. ABSCc) > ABSCd) 

Compute q = d/c 

thoi REAL (Zj./Za) = Ca 4^ b*q)/Cc ♦ d*q> 
IMAG CzjL/Za> = Cb - a*q)/{c ♦ d*q) 

2. ABSCc) < AESCd) 

Ca ♦ bI)/Cc ♦ dl) = Cb - aI)/Cd - cl) , 
which reduces to the first case. 

Hie comparison between ABSCc) and ABS Cd) 
is adequately performed in short preci- 
sion in both modules. 

Error and Exceptional Conditions: 

I : OVERFLOW, UNDERFLC^ and ZERODIVIDE 
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Positdiye Iii1:eqer Exponentia-tioni C fixed 
binary) 



Module Name: 



IHBWXItI 



Entry Point: IHEXIUO 

Function: 

To calculate z**©, where n is a ■ positive 
integer less than 2**31. 

Method : 

The contents of the target field are set 
to the value of z. The final result is 
obtained 1^ repeated squaring of the con- 
tents of the target field or squaring and 
multiplying hf x^ Multiplication is per- 
formed by the complex multiplication rou- 
tine IHEWfZU. 

Range: 

< n < 2**31. 

TSie precision rules of PL/I impose a 
furtter restriction in that if z has a 
precision Cp, q), this module may only be 
called if n*(p ♦ 1) - 1 < 31. This 
implies that n < 32/ Cp ♦ 1) < 16 for all 
such cases. 

Positive Integer Exponentiation (fixed 
decima 1) 



Module Name: 
Entry Point: 
Function: 



IHBWXIV 
IHEXIVO 



To calculate z**n, where n is a positive 
integer less than 2* ♦SI. 

Method: 

The contents of the target field are set 
to the value of the argument. The final 
result is obtained by repeated squaring 
of the contents of the target field or 
squaring and multiplying by the argument. 
Multiplication is performed b^ the con^- 
plesc multiplication routine IHEWMZV. 

Range: 

Ihe precision rules of PL/I impose the 
restriction that if z has a precision 
CPfq>f this BKJdule may only be called if 
n*Cp ♦ 1) - 1 < 15. This implies that 
n < 16/ Cp + 1) < 8 for all s\ich cases 
and, in fact, this module will operate 
only for the range < n < 8. 

Integer Exponentiation C floating-point) 

Module Names aiKi Entry Points: 



Module Entry 

Argument Name Point 

Stort float IHEWXIW IHEXIWO 

Long float IHEWXIZ IHEXIZO 

Functi on : 

To calculate z**nt where n is an integer 
between -2** 31 and 2 ♦♦SI - 1 inclusive. 

Method: 

If the exponent is and the argumont 
non-zero, the result 1 is returned imme- 
diately. If the exponent is non-zero, 
the contents of the target field are set 
to the argument value. The exponent is 
made positive and the argument raised to 
this positive power by repeated squaring 
of the contents of the target field or 
squaring and multiplying by the argument. 
Multiplication is performed by a branch 
to the complex imi It i plication subroutine. 
Then, if the exponent was negative, the 
reciprocal of the result is taken, other- 
wise it is left unchanged. 

Error and Exceptional Conditions: 

P : z = with n < 

I : OVERFLOW, UND121FLOW 

Since x*^C-m), where m is a positive 
integer, is evaluated as l/Cx^^m), 
the OVERFLOW condition may occur when 
m is large and the UNDERFLOW condi- 
tion when X is very small. 

H 2 OVERFLOW or UNDERFLOW in complex mul- 
tiplication routine (IHEWMZW or 
IHEWMZ2) 

General Floating-Point Expona:itiation 

Module Names and Entry Points: 



Arqum^it 
Short float 
Long float 

Function: 



Module 

Name 
IHEWXXW 
IHEWXXZ 



Entry 

Point 

IHEXXWO 

IHEXXZO 



To calculate Zj.^^Zaf where z^ and Za are 
complex numbers of the same precision. 

Method: 

Wh-en Zi = 0, the result is returned if i 
REMi(Za) > and IMAGCZa) = 0. 

Otherwise, ZjL**Za is computed as 

EXPCza^LOOCz^)) , 

with the proviso that if IMAGCzj.) = 
then LOG CABS (zj.) ) is calculated by a call 
to the real LOG routine, not to the com- 
plex LOG routine. 



Chapter 2. Arithmetic: Complex O'perations 11 



Error and Exceptional Conditions: 

P : Zj^ = with either REALCz^) < or 



O : a. REAL Cz2 ♦LOG (2j.n > 17*1.673: error 
caused in IHEWEXS or IHEWEXL 

b, IHEWXXW: 

ABS (IMAG(Za*LOGC2jL)H ^ 2**18*pi: 
error caused in SIN routine 
(IHEWSNS) 

IHEWXXZ: 

ABS(IMAG(Za*LOGCZi>)> > 2**50*pi: 

error caused in SIN routine 

CIHEWSNL> 



FUNCTIONS WITH REAL ARGUMENTS 
ADD (Fixed decimal) 
Module Name: IHEWADD 
Entry Point: IHEADDO 

Function: 

ADD(Xj.,Xa»p#q> where k^ and Xa are real 
fixed- point decimal numbers, and Cpfq) is 
the required precision of the result. 



Method: 

If both arguments are non- 
the module IHEWAPD is used 
one with the larger scale 
it the scale factor of the 
convert it to precision 31 
itients are added together, 
used to convert the sum to 
precision and to assign it 
field. 



zero, a call to 
to shift the 

factor to give 
other, and 
The argu- 

and IHEWAPD is 
the specified 
to the target 



If one of the arguments is zero, the 
other is treated as the sum above. 

Error and Exceptional Conditions : 

H : FIXEDOVERFLOW or SIZE may occur in 
IHEWAPD- 

MAX^MIN 

Module Names and Entry Points: 



Argument 

Fi xed bi na ry 



Fixed deciirial 



PL/ 1 Module Entry 

Function Name Point 

MAX IHEWMXB IHEMXBX 

MIN IHEMXBN 

MAX IHEWMXD IHEMXDX 

MIN IHEMXDN 



Short float 



Long float 



MAX IHEWMXS IHEMXSX 

MIN IHEMXSN 

MAX IHEI^XL IHEMXLX 

MIN IHEMXLN 



Function: 

To find the maximum or the minimum of a 
group of arithmetic values. 

All arguments must have the same base, 
scale and precision. 

Method: 

IHEWMXB, IHEWMXS,IHEWMXL: The value of 
the current maximum or minimum is set 
to the value of the first argument ; it 
is then compared algebraically with the 
next argument and replaced by it if 
appropriate. The process is repeated 
until a test on the argument list 
indicates that all source items have 
been processed, %rtien the current value 
is stored as the result. 

IHEWMXD: The address of the current 
maximum op minimum is set to the 
address of the first argument i this 
argument is then compared algebraically 
with the next argument, and the address 
of the latter replaces that of the for- 
mer if appropriate. The process is 
repeated until a test on the argument 
list indicates that all source items 
have been processed, when the result is 
moved into the target field. 



FUNCTIONS WITH COMPLEX Al^UMENTS 



ADD (Fixed decimal) 



Module Name: 
Entry Point: 
Function : 



IHEWADV 
IBEADVO 



ADDCZifZatPf g) where z± and Za are com- 
plex fixed-point decimal numbers, and 
Cp,q) is the required precision of the 
result. 

Method : 

The real parts of each argument are added 
and the sum is assigned to the target 
field by using the real fixed decimal ADD 
module CIHEWADD). The imaginary parts 
are treated similarly. 

Error and Exceptional Conditions: 

H : FIXEDOVERFLOW or SIZE may occur in 
IHiWAPD. 
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MULTIPLY (fixed binaryl 

Module Name: IHEWMPU , ^ 

Entry Point: IlEMPUO 

Function: 

M0LTIPLYCzj.,Z2iPf<3^ where Zj. and z^ are 
complex fixed- point binary numbers, and 
Cp, q) is the required precision of the 

result. 



Method: 

Let the arguments be Zj_ = a + bl and 

22 = c ♦ dl. 

Then REiyLCzj.*Za) = a*c - b*d 
IMA,G(Zj.*Za) = b*c + a*d 

The real and imaginary parts of the pro- 
duct are computed. These numbers are 
then shifted to giire them the required 
scale factorCq). 

The results of the shifts are t^ted for 
FIXEDOVERFLOW and truncated by left 
shifts. 

Error and Exceptional Conditions: 

I : FIXEDOVERFLCW 



DIVIDE (fixed binary) 
Module Name: IHEWDVU 



MULTIPLY (fixed decimal) 
Module Name: IHEWMPV 
Entry Point: IHEMPVO 
Function: 



M0LTIPLY(zj.iZ2fPiq) where Zj. and z^ are 
complex fixed-point decimal numbers, and 
(p,q) is the required precision of the 
result. 



Method: 



Let z± - a * bl and z^ 



c + dl, then: 



REAL(Za*Z2) = a*c - b*d. 
IMIlG(Zj,*Z2> = b*c ♦ a*d. 

Hie real and imaginary parts are calcu- 
lated and then each is assigned to the 
target with precision (p,q) by separate 
calls to the entry point IHEAPDA of the 
decimal shift and assign module IHEWAPD. 

Error and Exceptional Conditions: 

H : FIXEDOVERFLCM or SIZE in IHEWAPD. 



JEiitry Point: 
Fwiction : 



IHED¥UO 



DIVIDE (z-i,Z2tPfq) where Zj_ and Z2 ^^^ 
complex fixed- point binary numbers, and 
(p,q> is the required precision of the 

result. 



Method: 

Let Zj. = a ♦ bl, and Za = c ♦ dl, tben: 

REAL(zjl/Z2) = Ca*c + b*d)/(c**2 ♦ d**2) 
IMAG(za/Z2) = (b*c - a*d)/(c**2 ♦ d**2) 

The expressions a*c ♦ b*d, b*c - a*d, and 
c**2 ♦ d**2 are computed with a precision 
of 63- The denominator, c**2 ♦ d**2 is 
shifted to precision 31 by either a right 
or left shift. 

Two calls are then made to an incorpor- 
ated subroutine which accepts a numerator 
and shifts it so that it has two insigni- 
ficant leading digits. It then divides 
by c**2 + d**2 and shifts the quotient to 
the required scale factor (q) . 



Error and Exceptional Conditions: 
I : FIXEDOVERFLOW or ZERODIVIDE 

DIVIDE {fixed decimal) 



Module Name: 



IHEWDW 



Entry Point: IHEDWO 

Function: 

DIVIDE(zi.,22iPf<3) %*i.ere z^ and Za are 
complex fixed- point decimal numbers, and 
(p,q) is the required precision of the 

result. 

Method: 

Let Zi ^ a + bl, and Za = c ♦ dl, then 

JREAL(zj./Z2) = (a*c ♦ b*d)/(c**2 + d**2) 
IMAGCzx/Za) = (b*c - a*d>/(c**2 + d**2> 

The expressions a*c ♦ b*d, b*c - a*d, and 
c**2 + d**2 are computed. Leading zeros 
are remolded from the denominator (c**2 ♦ 

d**2) by truncation on the left and a 

left shift if necessary. If the 
denominator is still more than 15 digits 
long it is truncated on the right to 15 

digits. 
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T%io calls are then made to an incorpor- 
ated subroutine which accepts a numerator 
and shifts it to precision 31 with 2 
leading zeros by calling IHEWAPD (via 
entry point IHEJyPDB) • It then divides by 
c**2 ♦ d**2 and calls IHEWAPD (via entry 
point IHEAPDAI to assign the quotient to 
the target field with the required preci- 
sion (p,q). 

Error and Exceptional Conditions : 

I : ZERODIVIDE 



Entry Point: 
Functions 



IHEABVO 



H 



FIXEDOVERFLOW or SIZE in IHEWAPD 



ABS (fixed binary) 
Module Name: IHIHAKJ 
Entry Point: IHEABOO 

Function: 

To calculate ABS(z) = SQRT(x**2 ♦ y**2), 
where z = x ♦ yl . 

Method: 

If X = y, result is x*SQRT(2). 
Otherwise, 

let XI = MAXCABS(x),ABS(y>) 

Yl = MIN(ABS(x),ABS(y)). 

Then ABS(z) is computed as 

X1*SQRT(1 ♦ (Y1/X1)**2), 

where the fixed binary calculation of 
SQFT (g) for 1 < g < 2 is included within 
the module. 

The first approximation to the square 
root is taken as 

g/(l*g) ♦ (H-g)/il, 

with maximum relative error 1.8*2**-10. 
One Ne%rton-Raphson iteration gives maxi- 
mum relative error 1.6*2**-20, and suf- 
fices if XI < 2**(15-q) where q is the 

scale factor of z. 

Otherwise a second iteration is used, 
with theoretical maximum relative error 
of 1.3*2**-40. 

Error and Exceptional Conditions: 

I : FIXEDOVERFLOW 
ABS (fixed decimal) 
Module Name: IHEHABV 



To calculate ABS(2) = SQRT (x**2 + y**2) 
%*ere z = x ♦ yl. 



MetlK>d: 

X and y are converted to binary, with 
appropriate scaling if either exceeds 9 
significant decimal digits. 

Let XI be the maximum, and Yl the mini- 
mum, of the absolute values of the two 
binary numbers thus obtained. 

Then if XI = Yl = 0, result is 
returned. Otherwise, an approximation to 
ABS(2) is coropited as 

X1*SQRT(1 ♦ (Y1/X1)**2), 

where the fixed binary calculation of 
SQRT(gl for 1 < g < 2 is included within 
the module. 

The first approximation to the square 
root is taken in the form 

A * B*(l + g) - A/(l ♦ g) 

with maximum relative error 2.17*10**-1, 
and one Newton -Raphson iteration then 
gives a value with maximum relative error 
2.35*10**-8. 

Multiplication by XI produces a value for 
ABS(2) which is rounded and converted to 
decimal, and this suffices if it has not 
more than 7 significant decimal digits. 
Otherwise, this approximation is scaled 
if necessary and used in a final Newton- 
Raphson iteration for SQRT(x**2 ♦ y**2) 
in decimal, with theoretical maximum 
relative error 2.76*10**-16. 

Error and Exceptional Conditions: 

I s FIXEDOVERFIjOW 



ABS (floating-point) 

Module Names and Entry Points: 



Argument 



Module 
Name 



Entry 
Point 



Short float IHEWABW IHEABWO 
Long float IHEWABZ IHEABZO 

Function: 

To calculate ABS(z) = SQRT (x**2 ♦ y**2), 
where z = x ♦ yl . 

Method: 
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Let 2 = X ♦ yl. If X = y = 0, answer is 
0. 

Otherwise let z^ = MAX Ci^BS Cx) ,ABS Cy>> 

and Za = MIN CABS Cx> ,ABS Cy) ) . 

Then the answer is computed as 



ABSCz) = Zi*SQRTCl -§► (z^/Zjt) ♦♦2) . 
Accuracy: 

IHEWABW 

^ ^ .« — . . — . — -^ — ^ — . ^_— — .. ^ 

I Arguments | Relatiwe Error | 
I I ♦10**6 I 

^ — ^^^„ — ^^_^. — , — ^^ ^ — ^„^„^ 

I Range | Distribution! RMS | Maximum | 

I Full range! Exponential I 0.833 | 2.02 | 

I I radially, | | | 

I I uniform | | | 

I [round origin] | | 

L .-.. >X .- -L— — . X . J 

IHEWAB2 

^ ., — , .,_^ — ^ . ^ 

I Arguments | Relative Error | 
I I ♦10**15 I 

|. ^^^^^ . _„^.^. ^„^ .^ 

I Range | Distribution! RMS | Maximum | 

^ ^ ^^^^ ^ ^ 

I Full range I Exponential | 0.828 | 3.38 | 

I I radially, | [ | 

I I uniform | | | 

I I round origin | j | 



Error and Exceptional Conditions: 
I : OVERFLOW 
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CHAPTER 3: MATHEMATICAL BUXLT-IN FUMCTIONS 



The Library supports all floating point 
aritlimetic generic functions and has separ- 
ate modules for short and long precision 
real arguments. Additionally, the Library 
has separate modules for short and long 
precision complex arguments where these are 
admissible. 

Since the calling sequence generated in 
compiled code is the same as that required 
for passing the same arguments to a PL/ I 
procedure, it is permissible to pass the 
names of any of the float arithmetic gener- 
ic functions as arguments between proce- 
dures, according to the normal rules for 
entry names. 

Any restrictions on the admissibility of 
arguments are noted under the heading 
■Error and Exceptional Conditions. • 

Error and Exceptional Conditions ; These 
cover conditions which may result from the 
use of a routine; they are listed in four 
categories: 

P ~ Programmed conditions in the module 
concerned. Programmed tests are 
made where this is not too costly 
and, if an invalid argument is 
found, a branch is taken to the 
entry point IHEERRC of the execution 
error package CEXEP) . This results 
in the printing of an appropriate 
message and in the ERROR condition 
being raised- 

I — Interruption conditions in the 

module concerned. For those rou- 
tines %^ere an OVERFLOW interruption 
may occur, the condition is passed 
to the ON condition error handler 
CIHEWERR} and is treated in the 
normal way. For those routines 
where an UNDERFLOW may occur, the 
condition is disabled and both 
intermediate and terminal underflows 
are accepted as true zero. In cer- 
tain circumstances, however, %rtiere 
intermediate underflow may cause 
severe deterioration in the acciuracy 
of the result, the condition is 
avoided by programmed tests. 

O — Programmed conditions in modules 
called by the module concerned. 
These occur when invalid arguments 
are detected in the module called. 

H — As I, but the interruption condi- 
tions occur in the modules called by 
the module concerned. 



Accuracy 

In order to appreciate properly the 
ipeaning of the statistics for accuracy 
given with each module, some consideration 
of the limits and implications of these 
statistics is required. Because the size 
of a machine word is limited, small errors 
may be generated by mathematical routines. 
In an elaborate computation, slight inac- 
curacies can accumulate and become large 
errors. Thus, in interpreting final 
results, errors introduced during the 
various intermediate stages must be taken 
into account. 

The accuracy of an answer produced by a 
routine is influenced by two factors: (1) 
the accuracy of the argument and (2) the 
performance of the routine. 

Most arguments contain errors. An error 
in a given argument may have accumulated 
over several steps prior to the use of the 
routine. Even data fresh from input con- 
version may contain slight errors. Ihe 
effect of an argument error on the accuracy 
of an answer depends solely on the nature 
of the mathematical function involved and 
not on the particular coding by %riiich that 
function is computed within a routine. In 
order to assist users in assessing the 
accuHMlation of errors, a guide on the pro- 
pa gational effect of argument errors is 
provided for each function. Wherever pos- 
sible, this is expressed as a simple 
formula. 

The performance statistics supplied in 
this document are based upon the assumption 
that the arguments are perfect (i.e., 
without errors, and therefore having no 
argument error propagation effect upon an- 
swers) . Thus the only errors in answers 
are those introduced by the routines 
thonselves . 

For each routine, accuracy figures are 
given for the valid argument range or for 
represoitative segments of this. In each 
case the particular statistics given are 
those most meaningful to the function and 
range under consideration. 

For example, the maximum relative error 
and the root-mean- square of the relative 
error of a set of answers are generally 
useful and revealing statistics, but are 
useless for the range of a function where 
its value becomes 0, since the slightest 
error of the argument value can cause an 
unbounded fluctuation in the relative roag- 
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ni tilde of the answer. Such is the case 
with SIN(x) for values of x close to pi; in 
this rai^e it is more appropriate to dis- 
cuss absolute errors . 



The results were deriired from randCHii 
distributions of 5000 arguments per seg- 
ment, generated to be either uniform or 
exponential, as appropriate. It must be 
emphasized that each lvalue quoted for the 
maximum error refers to a particular test 
using the method described above, and 
should be treated only as a guide to the 
true maximum error. 

This explains, for example, why it is 
possible that the maximum error quoted for 
a segment may be greater than that found 
from a distribution of different arguments 
over a larger range which includes that 
segment. 



Hexadecimal Truncation Errors 

While the use of hexadecimal numbers in 
System/360 has led to increased efficiency 
and flexibility, the effect of the variable 
number of significant digits carried by the 
floating-point registers must be noted in 
making allowance for truncation errors. In 
the production of the PL/I Library, special 
care was taken to minimize such errors, 
whenever this could be accomplished at 
minor cost. As a result, the relative 
errors produced by some of the Library rou- 
tines may be considerably smaller than the 
relative error produced in some instances 
by a sir^le operation such as 
mult iplication . 

Representations of finite length entail 
truncation errors in any number system. 
With binary normalization, the effect of 
truncation is roughly uniform. With hexa- 
decimal normalization, however, the effect 
varies by a factor of 16 depending on the 
size of the mantissa i in a chain of compu- 
tations, the worst error committed in the 
chain usually prevails at the end. 

In short- precision representation, a 
nuii4>er has between 21 and 2'» significant 
binary digits. Therefore, the truncation 
errors range from 2**-2«* to 2**- 20 (5.96* 
10**-8 to 9.5*10**-7>. Assuming exact 
operands, a product or quotient is correct 
to the 24th binary digit of the mantissa. 
Hence truncation errors contritated by mul- 
tiplication or division are no more than 
2**-20. The same is true for the sum of 
two operands of the same sign. Subtrac- 
tion, on the other hand, is the commonest 
cause of loss of significant digits in any 
number system. For short -precis ion opera- 
tions, therefore, a guard digit is provided 
li^ich helps to reduce such loss. 



In long-precision representation, a 
number has between 53 and 56 significant 
binary digits. Therefore truncation errors 
range from 2**-56 to 2**-52 C1.39*10**-17 
to 2.22*10**-16). 

Normal care in numerical analysis should 
be exercised for addition and subtraction. 
In particular, when two algorithms are 
theoretically equivalent, it usually pays 
to choose the one which avoids subtraction 
between operands of similar size. 



Hexadecimal Constants 

Many of the modules described below dis- 
criminate between algorithms or test for 
errors by comparisons involving hexadecimal 
oonstantsi it must be realized that where 
decimal fractions are used in the descrip- 
tions the fractions are only quoted as con- 
venient approximations to the hexadecimal 
values actually employed. 

Algorithms 

The algorithms are the methods by which 
the mathematical functions are ccanputed. 
The presentation of each algorithm is 
divided into its major computational steps, 
with the formulas necessary for each step 
supplied. Some of the formulas are widely 
known; those that are not so %«idely known 
are derived from more common formulas. The 
process leading from the common formula to 
the computational formula is sketched in 
enough detail so that the deri lotion can be 
reconstructed by anyone who has an under- 
standing of college mathematics ai^ access 
to the common texts on numerical analysis.^ 

Many of the approximations were derived 
by the so-called "minima x" methods. The 
approximation sought by these metlKJds can 
be characterized as follows: given a func- 
tion f Cx) , an interval I, the form of the 
approximation (such as the rational form 
with specified degrees), and the type of 
error to be minimized Csuch as the relative 
error) , there is normally a unique approxi- 
mation to f Cx) whose maximum error over I 
is the smallest among all pc»sible approxi- 
mations of the given form. Details of the 
theory and the various methods of deriving 
such approxin^tions are left to the 
reference. 



^Any of the modem numerical texts may be 
used as a reference. One suc^ text is A. 
Ralston 's A First Course in Numerical Anal- 
ysis C McGraw-Hill Book Cc«npany, Inc. , New 
York, 1965). Background information for 
algorithms that use continued fractions may 
be found in H. S. Wall*s Analytic Theory 
of Continued Fractions (D. VanNostrand Co. 
Inc., Princeton, N.J., 1948>. 
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Termiiiolociy 

Ma.K±mvm and root-mean- s<juare naXues for 
the re la ti ire and Ciitiere necessarfl tlie 
ateolute errors are given for ^ch w>dule« 
These are defined thus : 



Let f ixl = the correct vaMe for a 
function. 
gfxl = tbe result obtained fron the 
aodule in question 

Then the absolute error of the result is 

AJKCf Cxi - gCxI), 
and the relati'^e error of the result is 

Em at in} - gCxJI/fC^I), 

Let the number of sample results obtained 
be N; then tlie root^-n^an-square of the 
absolute error is 

SQlTcH^CABSIf |JC|^) - g Cx !_)>♦♦ 2)/N) , 

and the root-mean- square of the relative 
error is 

SQRTfEiCABSflfCXil - gCXi) )/f f2Cj_)1**2>/N> • 

The Library mathematical modules are 
summarimed in Figures 5 and 6. 



I Function 

|SQRT 

|EXF 

|LOG,LOG2,LOG10 

I SIN, COS,SIND,COSD 

I TAN, TAND 

I AT AN, AT AND 

I SINE, COSH 

|TANH 

lATANH 

lEJRF, ERFC 



Real Arguments 
Short Float I Long Float 



IHEWSQS 
IHEWEXS 
IHEWLNS 
IHEWSNS 
IHEWTNS 
IHEWATS 
IHEWSHS 
IBEWTHS 
IHEWHTS 
IHEIiEFS 



IHEWSQL 

IHEWEXL 
IHEWLNL 
IHEWSNL 
IHEWTNL 
IHEWATL 
IHEWSHL 
IHEWTHL 
IHEWIffL 
IHEWEFL 



Figure 5. Mathematical Functions With 
Real Arguments 



I |.. ^ ^^ 

I Function | Short Float 

I SQRT I IHEWSQW 

I EXP I IBEWEXW 

I LOG I IHEWLNW 
|SIN,COS,SINH,COSH | IBEMSNW 

I TAN, TANH | IBEWTNW 

|ATAN, ATANH | IHEWATW 



Complex Argximents 
_ — ^-„_ , «. 

Long Float 



IHEWSQE 
IHEWEXZ 
IHEWLNZ 
IHEWSNZ 
IHEWTNZ 
IHEWATZ 



FUNCTIONS WITH REAL ARGOMENTS 
SORT Cshort floating-point real) 
Module Maii^s IBEWSC^^ 
Entry Point: IHBSQSO 
Function s 

To calculate the scpare root of «• 
Methods 

If X = 0, SQRTCx) = 0. Otherwise, let 
X = 16**(2*|> - q)*f, 

iRiliere p is an integer, q = or i, and 

1/16 < f < 1. Then 

SQRTCx) = 16**p*4*#-q*SQHTCf) 

The first approximation, y^, of SQRTCx) 

is obtained by the hyperbolic fit 

y^j = 16 ♦♦p»4**-q*C 1.681595- 
1, 288973/CO. 8%08065+f )> 

This approximation attains the minimax 
relative error. Tte maximum relative 
error is 2**-5,7i«8. 

Two Newton-Raphson iterations then yield 

yj. = Cyj^4-x/y^>y2 

ya = Cyj.-x/yi)/2+x/yi 

with a partial rounding* The maximum 
relative error of y^, is theoretically 
2**- 25.9. 

Effect of Argument Error: 

The relative error caused in the result 
is approximately half the relative error 
in the argument. 

Accuracy: 



Arguments | Relative Error j 

I ♦lo^^e I 

.„_^ ^ ^ ^ 

I I Range | Distribution! RMS | Maximum | 



I 



I Full Range I Exponential | 0.133 | 

i X . X ^-X- 



0.477 I 



Figure 6. Mathematical Functions With 
Complex Arguments 



Error and Exceptional conditions: 
P : X < 

SQRT Cloncr floating-point real) 
Module Name: IHEWSQL 
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Entry Point: IHESQLO 



Entrf Pointi IHEE3CS0 



Function: 

To calculate the square root of x« 

Method: 

If x = 0, ^RTCx) = 0. Otherviset let 
X = 16**f2*p - q)*ff where p is an inte- 
ger, q = or 1, and 1/16 < f < 1. Then 

SQBTCx) = 16**p*4»**'-q ♦SQRTff). 

The first approximation of SQRTCf) is 
computed as: 

y = 16**p*l»**(l~q)*0.2202Cf*0,2587} 

This approximation was chosen in order to 
permit the use of single precision 
instructions in the final iteration by 
making the quantity x/ya-ya below less 
than 16**Cp-8> . 

Four Newton -Raphs on iterations of the 

form y = Cyn ♦ x/y^) /2 are then ap- 

plied, two in short precision and two in 
long precision, the last being computed 

as 

SQKPCx) = ya ^- Cx/y3 - y3>/2 

with an appropriate truncation maneuver 
to obtain virtual roundii^. 

Hie maximum relative error of the final 

result is theoretically 2**-63.23. 



Effect of an Argui^nt Error: ■ 

The relative error caused in the result 
is approximately half of the relative 
error in the argum^it. 



Accuracy: 

^. , ^ — _-.-.-. — -._« — -.^ --« ^ ^^. — ^ 

I Arguments | Relative Error | 

I I *10**15 I 

I Range | Distribution | RMS | Maximtim | 

I Full range I Exponential | 0.0310| 0,109 I 

t ^ Jl . .x-.. X J 



Error and Exceptional Conditions; 
P : X < 

EXP (short floatinci- point real) 
Module Name: IHEWEXS 



Function: To calculate e to the power x. 

Itethals 

' If X < ^180.218, a zero result is 
returned immediately. 

Otherwise EXPfx) is calculated as 
follows : 

1. Divide X by LOG (2) and write 

y = x/.IjOG(2) = 4*a-b-d 

%*iere a and b are integers, < b< 3 and 
< d < 1. 

Th«i EXP(x) = 2**y = 16**a*2**-b*2**-d 

2. Conpute 2**-d by the following frac- 
tional approximation: 

2»*-d = l-2*d/C0.034657359*d**2+d-i- 
9, 951159*18-617. 97227/ Cd**2+87.ill7**97)) 

This formula can be obtained by the 
transformation of the Gaussian continued 

fraction 

EXPC-2)=l-2/Cl+z/C2-2/C3*2/C2-2/C5*2/C2-m 
/C7+Z/2-.. •))))>)> 

The maximum relative error of this ap- 
proximation is 2**-29. 

3. Multiply 2**-d by 2**-b 

<l. Finally multiply by 16**a by adding a 

to the characteristic of the result 
of step 3. 



Effect of Argument Error: 

The relative error caused in the result 
is approximately equal to the absolute 
error in the argument, i.e., to the argu- 
ment relative error multiplied by x. 
Bius for large values of x, even the 
round-off error of the argui^nt causes a 
substantial relative error in the answer. 



Accuracy: 



Arguments 



I Relative Error | 

I I ♦lO^^e 1 

I Range |Distribution| RMS | Maxijnuro | 

|-1 < X < 1| Uniform | 0.129 | 0.44i| | 



|Full Range! Uniform | 0.115 | 

L-- , .-J. — . . J «.-- .A- 



0.459 
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Error and Eacceptioiial Conditiors : 

I : OVERFLCH if M > l?t|.673 

EXP Clonq floa1:inc|-> point real) 
Module Name: IHEMEXL 

Entry Points IHEEXLO 

Function: To calculate e to the power 3c. 
Method: 



If X < -180.2187, return zero as the 
result. 

Otherwise EXPC-x) is calculated as 
follows: 

1. Divide x by UOGil} and let 

y=x/LOGC2) = 4*a-b-c/16 

%tiere a, b, and c are integers, < b < 
3, and < c < 15. Then, as an exact 

representa-tion for x, obtain 

X = C4*a-b-c/16)*LOGC2)-d 

%^ere the remainder d is in the range < 
d < BDG(2)/16. This reduction is carried 
out in extra precision. Then 

EXP(x> = 16**a*2**-b*2**C-c/16>*EXPC-d) 

2. Compute EXP{-d) by using a minima x 
polynomial approximation of degree 6 
over the range < d < LOGC2)/16. 
•Bie coefficients of this approxima- 
tion were obtained by taking the 
minimax of relative errors under the 
constraint that the constant term 
shall be exactly one. The relative 
error is less than 2**~56.87. 

3. Multiply EXPC-d) by 2**C-c/16), then 
halve the result b times. 

4. Finally, multiply by 16*^a by adding 
a to the characteristic of the result 
of step 3. 

Effecrt of an Argtmient Error: 

The relative error caused in the result 
is approximately equal to the absolute 
error in the argument, i.e., to the argu- 
ment relative error ma It i plied by x. 
Thus for large values of x, even the 
round-off error of the argument causes a 
substantial relative error in the answer. 



ilccuracy: 

j Arguments 



I Relative Error | 
I I ♦lO^^lS I 

I Range | Distribution! RMS | Maximum \ 

|-1 < X < 1| Uniform | 0.0543| 0.209 ( 

|„^ ^^^ „- +.« + i 

I Full range I Uniform | 0.0a72| 0.426 | 

I ^^X ^ X X J 

Error and Exceptional Conditions: 

I : OVERFLOW if x > 174.673 

LOG, LOG2, LQGIO (short floating-point 
real) 

Module Names IHEWLNS 

Entry Points: 

Matheaatical PL/I Entry 

Function Name Point 

Log X to the base e LOG(x) IHELNSE 

Log X to the base 2 LOG2 (x) IHELNS2 

Log X to the base 10 LOGlO(x) IHELNSD 



Function: To calculate log x. 

Method: 

Let X = 16**p*2** C-q)*m yhere p and q are 
integers, < q < 3, and 1/2 S m < 1. 

Two constants, a {= base point) and b (= 
-LOG2 Ca) ) , are defined as follows : 



If 1/2 < m < 1/SQRTC2); 
= 1 



then a = 1/2, b 



If 1/SQRT(2) < m < 1: then a = 1, b = 

Let y = Cm-a)/(m*a). 

Then m = a*Cl*y)/ Cl-y) and ABSCy) < 
0.1716. 

Now X = C2** C4*p-g-b))*CCl*y)/Cl-y)). 

•Eheref ore 

LOGCx) = C4*p-q-b)*LOGC2) ♦ 
LOG(Cl+y)/Cl-y)). 

To obtain LOG C (l+y)/Cl-y) ) first w = 2*y 
= (m-a)/(0.5ro*'0. 5a) is ccanputed Cwhicb is 
represented in System/360 with more sig- 
nificant digits than y itself), then the 
following approximation is performed: 

LOGC(l*y)/(l-y)) = w* (Cq 
'•^Ci*w«*2/ Cca-w**2) ) 

The coefficients were obtained h^ the 
minimax rational approximation of L0GCC1+ 
y)/Cl-y))/C2*y># ^^ relative error, under 
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the €X>ns1:raint that the first term shall 
be one- The maxipuin relative error of 
this approximation is less than 
2**-25.33. 



IiOG2Cx) or LOGlOCx) is calculated by mul- 
tiplying the above result by LOG2(e) or 
LOGlO(e) respectively. 



Effect of Argument Error: 

The absolute error caused in the result 
is approximately equal to the relative 
error in the argument. Thus if the argu- 
ment is close to 1, even the round-off 
error of the argument causes a substan- 
tial relative error in the answer, since 
the function value there is very small. 

Accuracy: 

^ , . ^^ , — ^^^ ^ . — ^ 

I Arguments | Relative Error | 
I I ♦10**6 I 

I — — . — ^ — .^^^ — _^. ^„_^ ^ 

I Range [Distribution! RMS | Maximum | 

i, . - X . X- .-X— . — . J 

IHELNSE 

J ^ , .^-.-. ^ ^ ^— I 

[Excluding | III 

|0.5 < x lExponential | 0.122 | 0.841 | 
l< 2.0 I III 

L X «.X ^ X ^ . J 

IHELNS2 

J. -.-^ , . ^ — ^ ^ . ^ 

I Excluding | II I 

|0.5 < x lExponential | 0.3*10 | 0.980 | 
|< 2.0 I II I 

L , X . -X X--- . S 

IHELNSD 

J. -^, — ^ , -^ ^ ^ , — —J 

I Excluding | III 

1 0.5 < X lExponential | 0.219 | 1.10 | 
l< 2.0 I II I 

L^ ^ . X .^ . -X , . X„- , J 

J. . ^ . , -^ — ,^. -J 

I Arguments | Absolute Error | 
I 1 ♦lO+^S I 

^ ^^^„„^„ ^^^ ^^^ 4 

I Range | Distribution] RMS | Maximum | 

L . X . _X, . X ^ J 

IHELNSE 

|0.5 < X I Oniform | 0.0255| 0.0679| 

|< 2.0 I II I 

L L X , X J 

IHELNS2 

1 0.5 < X I Uniform | 0.228 | 0.479 | 
l< 2.0 I III 

L X_ ^__X . X ^_J 

IHELNSD 

1 0.5 < X ( Uniform | 0.0228| 0.0720 1 
|< 2.0 I II I 

L X . ™X X-. J 



Error and Exceptional Conditions s 

P : X < 

IjOG> M)G2, LOGIO (long floating-point real) 

Module Name: IHEWLNL 

Entry Points: 

^fethematical PL/I Entry 

Function Name Point 

Log X to the base e LOGCx) IHELNLE 

Log X to the base 2 L0G2(x) 1HELNL2 

Log X to the base 10 LOGIO (x) IHELNLD 



Function: To calculate log x. 

Method: 

Let X = 16**p*2**(-q)*m where p and q are 
integers, < q < 3, and 1/2 < m < 1. 

Two constants t a C= base point) and b C= 
-LOG2(a)), are defined as follows: 

if 1/2 < m < 1/SQRT(2): then a = 1/2, 
b = 1 

if 1/SQRT(2) < m < 1: then a = 1, b = 0. 

1^^^ y = Cm - a)/Cm ••• a). 

Then m = a*(l ♦ y)/Cl - y) and ABS(y) 
< 0.1716. 

Now X = 2**C4*p - q - b)*Cl + y)/(l - y) 
Therefore 

LOGCx) = C«f*p - q - b)*LOGC2) + 
LOG C CI ♦ y)/Cl - y)). 

To obtain LOG C Cl+y)/Cl-y) ) first w = 2*y 
= Cm-a)/C0.5m+0. 5a) is computed C which is 
represented in System /360 with more sig- 
nificant digits than y itself), then the 
following approximation is performed: 

L0G((l4-y)/(l-y)) = w* (Co+c i*w**2/ 

{ca-w»*2)) 

The coefficients were obtained by the 
miniroax rational approximation of LOG (CI* 
y)/Cl-y))/ C2*y), in relative error, over 
the range y**2 < 0.0 2944 under the con- 
straint that the first term shall be 1. 
The maximum relative error of this ap- 
proximation is less than 2**-60.55. 

LOG2Cx) or LOGlOCx) is calculated by mul- 
tiplying the above result by L0G2Ce) or 
LOGIO ( e) respectively. 

Effect of an Argument Errors 

13ie absolute error caused in the result 
is approximately equal to the relative 



Chapter 3. Mathematical: Functions with Real Arguments 21 



error in the argtment. ^tma if the argu- 
ment: is close to 1, even the round-off 
error of the argument causes a substan- 
tial relative error in the answer, since 
the function value there is very small. 



Accuracy: 



Entry Points s 



Mathematical 


PL/I 


Entry 


Function 


Name 


Point 


SinCx radian si 


SIN Cx) 


IHESHSS 


SinCx degrees) 


SiNDCx) 


IHESl^Z 


CosCx radians) 


COS(x) 


IHESHSC 


Cos Cx degrees) 


COSD(x) 


IBESNSK 



^. , . _^_. — . . .» — . — . ^ 

I Arguments | Relative Error | 
1 I «10**15 I 

I , ^ — ^^„_^^„^. — — ^„„„«„4 

I Range | Distribution! RMS | Maximuro | 

L --.«_ L ^ . . «-.X ^-— . X—^ . J 

IHELNLE 

^^^, — ^. __^. .«-.^___ — ^-._^™.-.-_ ^ 

I Excluding | j j | 

|0.5 < x [Exponential j 0.05%4| 0.339 | 
l< 2.0 I II I 

L ._ . X- , . -X— .^ ^-X™. J 

IHELNL2 

j.„ . j.^. _. ^., _ — ^ — ^ — _ — _^ 

I Excluding | III 

1 0.5 < X (Exponential | 0.0881 | 0.425 | 
|< 2.0 I 11 I 

I , X , ^«._— — — X - X^- . J 

IHEI.NLD 

J. . J «., — .-.>_™«.^ ^ .^ . ^ 

I Excluding | ill 

I 0.5 < X lExponential | 0.0659 1 0.322 | 
l<__2.0 -^^1 i i 1 



Arguments 



I Absolute Error | 
I I ♦lO^^lS I 

I Range j Distribution! RMS j Maximum j 

L-__ X ^ . ^_X_ X_ . J 

IHELNLE 



|0.5 < X 
|< 2.0 

L , . 

IHELNL2 



I Unifomi 
I 



I 0.02391 

I I 

.X - X- 



0.0«I72| 

I 

J 



J. «._^ ^ . ^ . — ^^^ . ^ 

JO.5 < X I Uniform | 0.0291| 0.0576 | 
|< 2.0 I II I 

L_ . ^ X, ^ ^ - X— ^ X-« ^ J 

IHELNLD 



^ ^ , ^-. ^ — . ^ 

1 0.5 < X I Uniform 1 0.0125 | 0.0294 1 
|< 2.0 I II I 



Error and Exceptional Conditioi^ z 
P : X < 



sm, SIND, COS, COSD C short floatinqf- point 
real) 



Function: To calculate sin x or cos x* 
Method: 

Let k = pi/4 

Evaluate p = ABSfx)*(l/k) if x is in 

radians 
or p = ABSCx)* (1/45) if x is in 

degrees g 
using long-precision multiplication to 
safeguard accuracy. 

Separate p into integer part q and frac- 
tional part r, i.e., p = g ♦ r i4iere 
< r < 1. 



Define qj. 
qx 
qx 



q if SIN or SIND is required 

and X ^ 0| 

q ♦ 2 if COS or COSD is 

required! 

q + 4 if SIN or SIND is 

required and x < 0- 



Then for all values of x each case has 
been reduced to the computation of SIN(k* 
Iqi+r)) = SINCt) say, where t 2:: 0. 



Let 


qa 


= 


MOD 


Cqi,8). 






If 


qa 


= 


0, 


SIN(t) 


= 


SINCk^r) 


If 


qa 


= 


1, 


SINCt) 


= 


COSCk*Cl-r)) 


If 


qa 


= 


2, 


SIN(t) 


= 


COSCk*r) 


If 


q^ 


= 


3, 


SINCt) 


rr 


SINCk»Cl-r)) 


If 


qa 


=: 


^, 


SINCt) 


= 


-SlNCk*r) 


If 


^:i 


= 


5, 


SINCt) 


= 


-COSCk*ll-r)) 


If 


qa 


= 


6, 


SINCt) 


= 


-COSCk*r) 


If 


qa 


= 


7. 


SINCt) 


= 


-SINCk*Cl-r)) 



Module Names IHEHSNS 



Thus it is necessary to compute only SINC 
k^ra.) or C0SCk*rj.) where rj. = r or 1 - r 
and £ r^ < 1, as follows: 

1. SINCk*ri) = r3.*Cao 

♦aj.rj^**2^aarj.**4+a3rj.**6) 

The coefficients were obtained by the 
Chebyshev interpolation. The maximum 
relative error is less than 2**-28.1. 

2- cosCk^rj.) = l+bj,rj.*'»2*barj.**4*b3rj.**6 

The coefficients were obtained by a 
variation of the minimax approximation 
which provides partial rounding for 
the short precision computation. The 
maximum absolute error is 2**-24.57. 

Effect of an Argument Error: 
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The absolute error of the answer is 
approximately equal to tlie absolute error 
in the argument. Hence, the larger the 
argument, the larger its absolute error 
and the larger the absolute error of the 
result. Since the function diminishes 
periodically for both sine and cosine, no 
consistent control of the relative error 
can be maintained outside the range -pi/2 
to pi/2 radians (or -90 to +90 degrees). 

Accuracy: 



Arguments 



I Range 

L ^-_-. 



I Absolute Error 
I *10**6 

I Distribution I RWS | Maximum 

L , ^ — .„.-- -JL . J. 



IHESNSS 



^ — ^ -^ , -™^ ^ 

JABSCx) < I Uniform | 0.0467| 
I pi/2 I I I 

I pi/2 < I II 

|ABS(x) I Uniform | O.OilOO| 
1^ 10 I I I 

110 < i I I 

|ABSCx> I Uniform | 0.0«*01| 

|< 100 I II 

I -^ -L- . -«X-— ^ A. 

IHESNSC 

J. . ^ ^ ^ ^ ^. 

|0 < X < pij Uniform | 0.0408| 

1-10 < X I I I 

|< 0, I Uniform | 0.0il02| 

I pi < X I I I 

|< 10 I I 

|10 < I I I 

|ABS(x> I Uniform | 0.0398| 

1^ 100 I I I 

L _-.JL .^ .-., X-™^ JL. 

Error and Exceptional Conditions: 

P : IHESNSS, IHESNSC: 

ABSCx) > 2**18*pi 



0.119 i 
I 

-• i 

I 
0.125 I 

I 

-I 

I 
0.12*4 I 

I 



0.119 I 



0.120 I 
I 

I 

I 
0.113 I 

I 



Function: To calculate sin x or cos x. 



== ABSCx>/Cpi/4) for X in radians, 
- ABS(x)/45 for X in degrees, 
= q ♦ r, q integral, < r < 1. 



q for SIN or SIND with positive 

or zero argument, 

q ♦ 2 for COS or COSD, 

q -I- 4 for SIN or SIND with 

negative argument, 
and qa = MODCqi,8). 



Since COS(x) = SIN CABS Cx) -i- pi/2) 
and SINt~x> = SINCABSCx) + pi). 



Method: 




Let 

or 

and 


y 
y 
y 


Take 


q± 




qi. 



it is only necessary to find 



SINCpi/4*Cqa + r)), for < q^ < 7. 



Therefore compute; 



SINCpi/4*r), if qa = or 4, 

COSCpi/4*(l - r)}, if qa = 1 or 5, 

COSCpi/4*r), if qa = 2 or 6, 

SIN{pi/i|*Cl - r)) if qa = 3 or 7. 



SINCpi/4*rj.)/ri, where r^ is r or 
CI - r), is computed by using the Cheby- 
shev interpolation polynomial of degree 6 
in ri**2, in the range < rj^**2 < 1, 
with maximum relative error 2**C-58). 



C0S(pi/1*rj.) is computed by using the 
Chebyshev interpolation polynomial of 
degree 7 in rj.**2, in the range < r^.** 
2 < 1, with maximum relative error 

2**C-6«*.3). 



IHESNSZ, IHESNSK: 
ABSCx) > 2**18*180 



Finally, if qa ^ ^ a negative sign is 

given to the result. 



SIN, SIND, cos, CQSD Clonq floating- point 
real> 



Effect of an Argiiment Error: 



Module Name: 



Entry Points : 



IHEWSNL 



Mathemati cal 

Function 
SinCx radians) 
SinCx degrees) 
CosCx radians) 
CosCx degrees) 



The absolute error of the answer is 
approximately equal to the absolute error 
in the argument. Hence, the larger the 
argument, the larger its absolute error 
PL/I Entry and the larger the absolute error of the 
Name Point result. Since the function diminishes 
SINCx) IHESNLS periodically for both sine and cosine, no 
SINDCx) IHESNLZ consistent control of the relative error 
COSCx) IHESNLC can be maintained outside the range -pi/2 
COSDCx) IHESNLK to pi/2 radians Cor -90 to +90 degrees). 
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Accuracy: 

IHE3NLS 

I Arguments 



I Relative Error 
I ♦10**15 



I Range 

i . . 

IHESNLS 



I Distribution! RMS | Maximum 

.J. ^«x x^. . -.-. 



j ABSCx) < I j t i 

I pi/2 I Uniform | 0.0181 I 0.0771 | 
|. ^ + 4.^„„ „+ ^ 

I pi/2 < I II I 

I ABSfx) I Uniform j 0.317 | 2.36 | 

LL^^ fill 

I 10 < I I I I 

I ABSCx) I Uniform | 0.928 | 2.65 | 
I < 100 I I I I 

i . _X . ^-JL . A . . J 

IHESNIiC 

I < X < pi I Uniform | 0.0739 | 0.266 | 

I -10 < X I I I I 

I < 0, I Uniform | 0.0683 | 0.266 | 

I pi < X I I I I 

I < 10 I I I I 

I 10 < I I I I 

I ABSCx) I Uniform | 1.02 | 2.68 | 
I < 100 I I I I 

L X„-^ X J. 1 



Error and Exceptional Conditions: 

P : IHESNLS, IHESNLCs 

ABSCx) > 2**50*pi 

IHESNLZ, IHESNLK: 
ABSCx) > 2**50*180 



Let q and r fce respectively the integral 
and fractional parts of p. 



If q is even, put s = r; 
if q is odd, put s = 1-r, 



Let q^. = M0DCq,4). Then 



If 
If 
If 
If 



qj- 



0, 

It 

2, 
3. 



TAN CABS Cx)) 
TAN CABS Cx)) 
TAN C ABSCx)) 
TAN C ABSCx)) 



TANCpi*S/i*) 

COTCpi*S/i|) 

-COTCpi*S/i|) 

-TANCpi*s/i|) 



Compute TANCpi*s/i|) and CX>TCpi*s/4) as 
the ratio of two polynomials: 



TANCpi*s/4) = s*PCu)/QCu) 
COTCpi*s/4) = QCu)/Cs*PCul> 



where u = s**2/2 and 



PCul = -8.160901+u and 



QCul = 10. 7727511+5. 703366*u 
-0.159321*u*«2 



These coefficients were obtained by the 
minimax rational approximation in rela- 
tive error of the above form. The maxi- 
mum relative error of this approximation 
is 2**-26. The variable u, rather than 
s**2, was chosen for P and Q in order to 
improve the rounding effect of the 
coefficients . 



TAN, TAND Cshort floating-point real> 
Module Name: IHEWTNS 

Entry Points: 



Ma theroat i cal 
Function 



PL/ 1 


Entry 


Name 


Point 


TANCx) 


IHETNSR 


TANDCx) 


IHETNSD 



TanCx radians) 
TanCx degrees) 



Function: To calculate tan x. 

Method: 

Evaluate p = Cil /pi) ♦ABSCx) if x is in 
radians , 
or p = C 1/45) ♦ABSCx) if X is in 
degrees, 

using long- precis ion multiplication to 
safeguard accuracy. 



Finally, if x < 0, put 



TANCx) = -TAN CABS Cx)) 



Effect of an Argument Error: 

The absolute error of the answer is 
approximately equal to the absolute error 
of the argument multiplied i^ CI "•■ TANCx) 
♦♦2) . Hence if x is near an odd multiple 
of pi/ 2, an argument error will produce a 
large absolute error in the answer. 



The relative error in the result is 
approximately equal to twice the absolute 
error in the argument divided by SINC2^ 
x) . Hence, if x is near a multiple of 
pi/2, an argument error will produce a 
large relative error in the result. 
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Accuracy: 



Arguments 



I Relative Error 
I ♦lO+^S 



I Range (Distribution j RMS j Maxiroum 

L ^. X----^ ._.„-.«X_.«__-.^X 



IHETNSR 



JABSCx) < 


T- 


Uniform 


1 0.290 1 1.64 


1 pi/14 

H— 

|pi/i| < 
lABSCx) 


-+- 


Uniform 


i 1 

t 1 

1 0.369 1 1.54 


1 < pi/2 

I pi/2 < 
lABSCx) 


-+- 




1 I 

1 1 

1 0.321 1 4.81 


Uniform 


|< 10 






1 t 


I- — 

110 < 


-+- 


. ... — 


1 1 


|ABS(x) 




Uniform 


1 0.310 1 1.38 


|< 100 






1 1 


L ^ , 


_X_ 




_ X X 






Error and Exceptional Conditions: 

P : IHETIiSR: ABSCx) > 2**18*pi 
IHETNSD: ABSCx) > 2**18*180 

I : IHETNSR: OVERFLOW 
IHETNSD: OVERFLOW 



compute TAIiCpi*s/4) and COTCpi*s/4) as 
the ratio of two polynomials: 



TAMCpi*s/4) = s*PCs**2)/QCs**2)) 
COTCpi*s/4) = QCs**2)/(s*PCs**2)) 



where both P and Q are polynomials of 
degree 3 in s**2. The coefficients of P 
and Q were obtained by the miniroax 
rational approximation Cin relative 
error) of TAN(pi*s/4) of the indicated 
form. The maximum relative error of this 
approximation is 2**-55.6. 

Finally, if x < , TAN Cx) = -TAN CABS Cx) ) . 



Effect of an Argument Error: 

The absolute error in the result is 
approximately equal to the absolute error 
in the argument multiplied by Cl-»-TANCx)** 
2). Hence, if x is near an odd multiple 
of pi/2, an argument error will produce a 
large absolute error in the result. 

The relative error in the result is 
approximately equal to twice the absolute 
error in the argument divided by SINC2* 
x) . Hence, if x is near a multiple of 
pi/2, an argument error will produce a 
large relative error in the result. 



TAN, TAND Clonq f loa ting -^ point real) 

Module Name: IHEWTNL 

Entry Poi nts : 

Mathematical PL/I Entry 

Function Name Point 

TanCx radians) TANCx) IHETNLR 

TanCx degrees) TANDCx) IHETNLD 

Function: To calculate tan x. 

I^ethod: 

Evaluate 

p = C4/pi) ♦ABSCx) if x is in radians 
or p = Cl/45)*ABSCx) if x is in degrees. 

Let q and r be respectively the integral 
and fractional parts of p. 

If q is even, put s = rj 
If q is odd, put s = 1 - r. 

Let q^ = M0DCq,4). Then 



If 
If 

If 
If 



qi 



0, TAN C ABSCx)) 

1, TAN CABS Cx)) 



= TANCpi*s/4) 
= C(3TCpi*s/4) 



2, TANCABSCx)) = -COTCpi*s/4) 

3, TANCABSCx)) = -TANCpi*s/4) 



Accuracy: 
IHETNLR 



j Arguments 

I 

I Range [Distribution 

lABSCx) < I Uniform 
I pi/ 4 

I pi/4 < 

lABSCx) I ♦Uniform 

|< pi/1.5 



h 

I pi/1. 5 < 

lABSCx) 

1< 10 

j. .^^ ^ 

|10 < 
t ABSCx) 

i< 100 






♦Uniform 






♦Uniform 



Relative Error 
♦10^*15 



^^ ^^ ^^^ ^ 

RMS I Maximum 



0.646 



0.471 



84.2 



78.8 



0.571 



2.26 



H 



4730 
2710 



I ♦The errors quoted are those encountered 
I in a sample of 5000 points? each figure 
I depends very much on the particular 
I points encountered near the singulari- 
I ties of the function, where no error 
I control can be maintained. 

L .- ^-^ . . ., ^ .. . .^ 



H 
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Error and Exceptional Conditions s 

P : IHETNLR2 ABS(x) > 2**50*pi 
IHETNIiD: ABS(x> > 2**50*180 

I : IHETMLR: OVERFLOW 
IHETNID: OVERFLOW 



(siK>rt floating-point real) 



Module Name: 



IHEWATS 



For the basic range ABS(x) < TAN (pi/12), 
use an approximation forronla of the form 



ATANCx)/x = a + b*x**2 + c/(d * x**2) 

with relative error less than 2**-27.1. 

3. ATANDCx) and ATAND(y,x) 

The treatment is as above with the addi- 
tion of a final conversion of the result 
to degrees. 



Entr^ Points: 

Mathanatical PL/I Entry 

Function Name Point 

Arctan x Cradians) ATANCx) IHEATSl 

ArctanCy/x) Cradians) ATANCy,x) IHEATS2 

Arctan x Cdegrees) ATAND(x) IHEATS3 

ArctanCy/x) Cdegrees) ATANDCy^x) IHEATSi* 



Function: 

To calculate arctan x or arctan Cy/x). 
The result range is: 

Arctan x Cradians) ± pi/2 
Arctan Cy/x) Cradians)! pi 
Arctan x Cdegrees) ± 90^ 
Arctan Cy/x) (degrees)! 180° 

Method: 

1. ATANCy,x) 

If X = or ABS(y/x) > 2**2«4 , the answer 
SIGN(y)*pi/2 is returned except for the 
error case x = y = 0. Otherwise 

ATAN(y,x) = ATAN(y/x) if x > 
or ATAN(y,x) = ATAN(y/x) + SIGN(y)*pi 
if X < 0. 

Hence the computation is now reduced to 
the single argument case. 

2. ATAN(x) 

The general case may be reduced to the 
range < x < 1 since 

ATANC-x) = -ATANCx), and 
ATANCl/ABSCx)) = pi/2 - ATAN C ABS C x) ) . 

A further reduction to the range ABSCx) < 
TAN C pi/1 2) is made by using 

ATANCx) = pi/6 + ATAN C CSQRTC3)*x - 1)/ 
Cx ^ SQRTC3) )). 

Care is taken to avoid the loss of signi- 
ficant digits in computing 

SQRTC3)*x - 1. 



Effect of an Argument Error: 

Let t = X or y/x; then the absolute error 
of the answer approximates to the abso- 
lute error in t divided by CI ♦ t**2) . 
Hence, for small values of t, the two 
errors are approximately the samei howev- 
er, as t becomes larger the effect of the 
argument error on the answer error 
diminishes. 

Accuracy: 

J. — . ^ — . ^ . ^ 

I Arguments | Relative Error | 
I I ♦lO^^e I 
^ ^ + ^ 4 

I Range | Distribution! RMS | Maximum | 

i X -X ^ X 1 

IHEATSl 

J. ^ ^ J. ^ , — «^ 

I ABSCx) I Uniform | 0.127 | 0.898 | 

l< 1 I II I 

^-- ^ ^-^^^ ^ 4 

I Full range) Exponential [ 0.2«l6 | 0.99i* | 

L , X ^ X X J 

IHEATS2 

P ^ ^ . ^ ^ 

lABSCy) I Exponential | | j 

|< 1, I- i 0.291 I 1.62 I 

I ABSCx) I Uniform | | | 

1^ 1 I II I 

i ^ X ^-X X J 



Error and Exceptional Conditions: 
P : IHEATS2, IBEATS**: x = y = 



ATANCX), ATANDCX), ATAN CY.X)^ ATAND CY^X) 
Clonq flcatinq- point real) 



Module Name: 



Entry Points: 



IHEWATL 



Mathematical 
Function 
Arctan x Cradians) 
Arctan Cy/x) Cradians) 
Arctan x Cdegrees) 



PL/I 

Name 
ATANCx) 
ATANCy,x) 
ATANDCx) 



Entry 

Point 

IHEATLl 

IHEATL2 

IHEATL3 



Arctan Cy/x) (degrees) ATAND CY,X) IHEATL4 
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Function: 



ATANDCx) and ATANDCy,x> 



To calculate arctan m or arctanCy/x). 
The result range is: 

Arctan k Cradians) ± pi/2 

Arctan Cy/x) (radians)! pi 

Arctan x (degrees) ± 90<' 

Arctan(y/x) (degrees)! 180° 



Method: 

1. ATAN(y,x) 

If X = or ABS(y/x) > 2**56, the answer 
SIGN(y)*pi/2 is returned except for the 
error case x = y = 0. Otherwise 

ATAN(y,x) = ATAN(y/x) if x > 
or ATAN(y,x) = ATAN(y/x) + SlGN(y)*pi 
if X < 0. 

Hence the computation is now reduced to 

the single argument case. 

2. ATAN(x) 

The general case may be reduced to the 
range < x < 1 since 

ATAN(-x) = - ATAN(x), and 
ATAN(l/ABS(x)) = pi/2 - ATANCABS (x) ) . 

A further reduction to the range ABS(x) 

< TAN (pi /I 2) is roade by using 

ATANCx) = pi/6 + ATAN ( (SQRT (3) ♦x - 1)/ 
(x + SQRT(3))) 

Care is taken to avoid the loss of signi- 
ficant digits in computing 

SQRT(3)*x - 1 

For the basic range ABS(x) < TAN(pi/12), 
use a continued fraction of the form 

ATANCx) /x = 1 + u* (bo~ai/(tj. + u-a2/(b2+u 
-aa/Cba+u)))) 



where u 



x**2. 



The relative error of this approximation 
is less than 2**-60.7. 

The coefficients of this formula were 
derived by transforming a miniroax ration- 
al approximation in relative error over 
the range < u < 0.071797 for ATAN(x)/x 
of the following form: 

ATAN(x)/x = ao+u*( CCo+Cx*u-i"Ca*u*u+ 

G3*u*u*u)/ (do-i-di*u+d2*u*u-»' 
u*u*u) ) • 



The treatment is as above with the addi- 
tion of a final conversion of the result 

to degrees. 

Effect of an Argument Error: 

Let t = X or y/xj then the absolute error 
of the answer approximates to the abso- 
lute error in t divided by CI + t**2) • 
Hence, for small values of t^ the two 
errors are approximately the same; howev- 
er, as t becomes larger the effect of the 
argument error on the answer errcr 
diminishes. 

Accuracy: 

^ — .-..- — . ^-. — «— .. — -.^ «., — .^^ — --^ .^ 

I Arguments | Relative Error | 
I I ♦lO^^lS I 

^^^„„^, — ^„ ^ — ^ ^^ ^ 

I Range | Distribution! RMS | Maximum | 

i -- . X , J. -X J 

IHEATLl 

lABSCx) I Uniform | 0.0415 | 0.206 | 
1< 1 t I I I 

^„^^ ^^^^ +„ ^„4 4 

JFull range I Exponential I 0.0526 | 0.206 j 

L— -. X-.. . X ^ X ^--J 

IHEATL2 



under the constraint that a 



1. 



lABSCy) I Exponential I | [ 

1^ It 1 I I I 

I |. — «^-. ^ 0.0688 I 0.358 | 

JABSCx) I Uniform | | j 

1^ 1 I I i I 

I , ^_— .X ^ . X ^-X I 

Error and Exceptional Conditions; 
P : IHEATL2, IHEATL4 : x = y = 

SINH, COSH (short floating-point real) 

Module Name: IHEWSHS 

Entry Points: 

Mathematical PL/I Entry 

Function Name Point 

Hyperbolic sin x SINHCx) IHESHSS 

Hyperbolic cos x COSH(x) IHESHSC 

Functions 

To calculate hyperbolic sin x or hyper- 
bolic cos X. 

Method : 

1. ABS(x) < 1 

Compute SINHCx) as: 

SINHCx)=x+Cj.*x**3+C2*x**5+C3*x**7 
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The ooef ficients were obtained by the 
miniiDax appxoxiroation Cin relative error) 
of SINHCxI/x as a function of x**2« The 
maxiiiiim relative error of this approxiDia- 
tion is 2** (-25.6) . 



Thus, for large values of x, even the 
round-off error of the argument causes a 
substantial relative error in the answer, 

Accuracy : 



x>l 



Compute SINH(x) as: 

SINH (x) = (1+D) * (EXP ( x*LOG (V) ) - 
V**2/EXP(X'i-L0G(V))) 

Using module IHEWEXS. 

Here l+jD=i/C2*V) , so that this expression 
is theoretically equivalent to (EXP(x)- 
EXF(-x))/2. The value of V (and conse- 
quently those of LOG(V) and D) was so 
chosen as to satisfy the following 
conditions : 

a) V is slightly less than 1/2 , so 
that D is positive and small 

b) IX>G(V) is an exact multiple of 
2#*(-i6). 

Condition (b) ensures that the addition 

x4-LOG(V) is carried out exactly. 

3. x<-l 

Use the identity 

SINH (x) =-SINH (ABS (x) ) 

to reduce to case (2), above. 

«♦. COSH(x) 

For all legal values of arguments, use 
the identity 

COSH (x) = (1+D) ♦ (EXP (x+LOG (V) ) ♦ 
V**2/EXP (X+LOG (V) ) ) 

Here the notation and considerations are 
identical to those used in the computa- 
tion of SINH(x), in (2) above. 

Effect of Argument Error: 

The relative error caused in the result 
is approximately as follows: 

SINH: The absolute error in the argument 
divided by TANH(x), i.e., of the 
order of the absolute error in the 
argument for large x, or of the 
relative error in the argument for 
smal 1 X « 

COSH: The absolute error in the argument 
multiplied by TANB(x), i.e., of the 
order of the absolute error in the 
argument. 



Arguments 



-^ ^ 

I Felative Error | 
I ♦lO^^e I 



I Range 

i ^ 



I Distribution I RMS | Maximum 

.x^ . i. X 



H 



IHESHSS 

J.-. J. 1 ^ .^^ 

|0 < I II I 

|ABS(x) < 1| Uniform \ 0.198 | 0.877 | 

|1 < I II I 

|ABS(x) < 2 1 Uniform | 0.255 | 1.03 I 

^ « + .^4. + — ^4 

I |ABS(x)<170 1 Uniform I 0.201 | 0.816 | 

I , X X X ^ I 

IHESHSC 

J. . J. — , J. ^ . ^ 

|ABS(x) < ij Uniform j 0.406 | 0.962 j 

|1 < I II I 

|ABS(x) < 2| Uniform | 0.248 | 0.720 | 
^ ^«-+ + 4^^ ^„4 

I |ABS(x)<170 I Uniform | 0.202 | 0.816 | 

L ^ X X X I 

Error and Exceptional Conditions: 

H : OVERFLOW in real EXP routine 
(IHEWEXS). 



COSH, SINH (long floating-point real) 

Module Name: IHEWSHL 

Entry Points: 

Mathematical PL/I Entry 

Function Name Point 

Hyperbolic cos x COSH(x) IHESHLC 

Hyperbolic sin x SINH(x) IHESHLS 

Functions 

To calculate hyperbolic sin x or hyper- 
bolic cos X. 

Method : 

1. ABS(x)<0. 8813711 

Compute SINH(x) as 

SINH(x) = Co*x+Cj.*x**3+Ca*x**5* 

4^C6*X**13 

The coefficients were obtained by the 
rainimax approximation (in relative error) 
of SINH(x)/x as a function of x**2. The 
maximum relative error of this approxima- 
tion is 2**(-55.7). 
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X20. 881374 



Compute SINH(x) as 

SINeCx)=Cl+D)*CEXPCx+LOG(V))- 
V**2/EXPCx'i"LOGCV)>) 

using module IHEEXL 

Here 1+D=1/C2*V) so that this expression 
is theoretically equivalent to (EXP(x)- 
EXPC-x))/2. The value of V Cand conse- 
quently those of LOGCV) and D) was so 
chosen as to satisfy the following 
conditions: 

a) V is slightly less than 1/2 so that 
D is positive and small. 

b) LOGCV) is an exact multiple of 
2**C-16). 

Condition Cb) ensures that the addition 
x*IX>G(V) is carried out exactly. 

3. x<-0.88137«l 
Use the identity 

SINK Cx) =-SINH CABS Cx) ) 
to reduce the case to that of step C2). 

4. COSHCx) 

For all legal values of arguments , use 
the identity: 

COSH Cx) = Cl + D) ♦ CEXP f X+LOG CV) ) * 
V**2/EXPCx*L0GCV))) 

Here the notation and considerations are 
identical to those used in the computa- 
tion of SINH Cx) in step C2) above. 



Effect of an Argument Error: 

The relative error caused in the result 
is approximately as follows: 

SINH: The absolute error in the argument 
divided by TANHCx), i.e., of the 
order of the absolute error in the 
argument for large x, or of the 
relative error in the argument for 
small X. 

COSH: The absolute error in the argun^nt 
multiplied by TANHCx), i.e., of the 
order of the absolute error in the 
argument. 



Thus, for large values of x, even the 
round-off error of the argument causes a 
substantial relative error in the answer. 



Accuracy: 

^-«- — --.-«.« , .^^^^ . ^ 

I Arguments | Relative Error | 
I I ♦lO^^lS I 
^^„„^ ^ , 4.^.„„„^„»^ ^ 

I Range | Distribution! RMS | Maximum | 

I .-. . A ^ . A ^ JL . ^ 1 

IHESHLC 

^ -..- ^ ^— , ^ . — ^ — , 

JABSCx) < 5 1 Uniform | 0.106 | 0.376 t 

^^^ ^^ ^ „ + + ^i 

JABSCx) I Uniform 1 | I 

|< 170 I I 0.109 I 0.390 I 

i-- ^ . X-. ^^ J.-^ X . ^ 1 

IHESHLS 

^^ .^ — .^^^y. ^ . — ^ — ^ ^— ^ ^ 

jABSCx) < I Uniform | 0.0373 | 0-203 | 
|0. 8813711 I III 

10.8813711 <| III 

lABSCx) < 5 1 Uniform | 0.100 | 0.354 | 

lABSCx) I Uniform | | | 

|S 170 I 1 0.102 I 0.361 I 



Error and Exceptional Conditions: 

H 2 OVERFLOW in real EXP routine 
CIHEWEXL). 



TANK C short floating-point real) 

Module Name: IHEWTHS 

Entry Point: IHETHSO 

Function: To calculate hyperbolic tan x. 

Method : 

1. ABSCx) < 2**-12 

Return x as result. 

2. 2**-12 < ABSCx) < 0.7 

Use a fractional approximation of the 
form: 

TANHlx)/x = 1 - 
x**2*C.0037828+.8i«l5651/Cx**2+2.i|717il9>) 

The coefficients of this approximation 
were obtained by taking the fdnimax of 
relative error, over the range x**2 < 
0.49, of approximations of this form 
under the constraint that the first term 
shall be !• The maximum relative error 
of this approximation is 2**-26.4. 

3. 0.7 S n < 9.011 

Use TmNHfxl = 1 - 2/lEXPC2*k) ♦ II. 
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i*. X ^ 9.011 

Return result 1. 

5. X < ~ 0.7 

Use the identity: 

TJy^HCx) = -TANBC-x)- 

and apply 3 or ** above t as appropriate. 

Effect of an Argument Error: 

The relative error caused in the result 
is approximately twice the absolute error 
in the argument divided by SINH(2*x). 
Thus for small values of x it is of the 
order of the relative error in the argu- 
ment, and as X increases the effect of 
the argument error is diminished. 

Accuracy: 

J, , ^ ^ , ^ 

I Arguments | Relative Error | 
I I *10**6 I 
^^ ^„^ ^ ^ 4 

I Range | Distribution! RMS | Maximum | 

|ABSCx><0.7 1 Uniform | 0.149 | 0.781 | 

|0.7<ABS(x) lUniform | 0.0389| 0.288 | 

|<9.011 I III 

L . X ^ ^±. X J 



«l. X ^ 20.101 

Return result 1. 

5. X < -0.5tl931 

Effect of an Argument Error: 

The relative error caused in the result 
is approximately t%iice the absolute error 
in the argument divided by SINiiC2«x). 
Thus for small values of x it is of the 
order of the relative error in the argu- 
ment, and as x increases the effect of 
the argument error is diminished. 

Accuracy : 



h 



. . ^ . . ^ 

Arguments | Relative Error 
I ♦10**15 



H 



I Range | Distribution! RMS | Maximum 

|ABSCx> < I Uniform | 0.0385| 0.192 
10.54931 I 1 I 



H 



10.54931 < I f I 

(ABSCx) < 5 1 Uniform | 0.0109| 0.160 

I ^-X X- X J 



ATANH Cshort floating-point real) 
Module Name: IHEWHTS 



TANH (long floating-point real) 

Module Name: IHEWTHL 

Entry Point: IHETHI.0 

Function: To calculate hyperbolic tan x. 

Method : 

1. ABS(x) < 2**-28 
Return x as result 

2. 2**-12 < ABSCx) < 0.54931 

Use a transformed miniroax approximation 
of the form 

TANH (x) /x=c o+d3L*x**2/ (x**2-i'Cj.) +da/ Cx«* 
2+C2)+d3/(x**2+C3) 

The roinimax of relative error was taken 
over the range x**2 < 0.30174 under the 
constraint that the first term is 1. 

The maximum relative error is 2**-63 

3. 0.54931 < X < 20.101 
TANH(x) = l-2/CEXPC2*x)+l) 



Entry Point: IHEHTSO 

Function: To calculate hyperbolic arctan x. 

Method: 

1. ABSCx) < 0.2 

Use a rational approximation of the form: 
ATANHCx) = X + X ♦♦ 3/ (a + b*x**2) 

2. 0.2 < ABS(X) < 1 

ATANHCx) = -SlGNCx)*0.5*LOGCC0.5 - 

ABSCx/2))/C0.5 ♦ ABSCx/2))) 



Effect of an Argument Error: 

The absolute error caused in the result 
is approximately equal to the absolute 
error in the argument divided by CI - x** 
2). Thus as X approaches +1 or -1, rela- 
tive error increases rapidly. Near 
X = 0, the relative error in the result 
is of the order of that in the argument. 
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Accuracy: 

1 Arguments 

I 

h 



I Relative Error 



I Range j Distribution | RMS | Maximum | 

|-0,2 < X I Uniform | 0-456 | 1.07 | 

|< 0.2 I III 

i« — . ^^ ^^^-.^-„^.^ — „^^ — ^ — »^ 

1-0.9 < X I Uniform | 0.391 | 1.18 | 
i< 0.9 I II I 

L , ^ .--.-.JL - .^^--JL™.- JL — ^ ^ ^J 

Error and Exceptional Conditions : 

P : ABSCx) > 1 

ATANH (long floating-point real) 
Module Name: IH^HTL 
Entry Point: IHEHTLO 

Function: To calculate hyperbolic arctan x. 
Method: 

1. ABS(x) < 0.25 



Use a Chebyshev polynomial of degree 8 in 
x**2 to compute ATANH(x>/x. 

2. 0. 25 < ABSCx) < 1 

ATANH(x) = -SIGN(x)*0.5*LOG'( C0.5 - 

ABSCx/2))/(0.5 ♦ ABS(x/2))) 

Effect of an Argument Error: 

The absolute error caused in the result 
is approximately egual to the absolute 
error in the argument divided by ( 1 - x** 
2). Thus as X approaches +1 or -1, rela- 
tive error increases rapidly. Near 
X = 0, the relative error in the result 
is of the order of that in the argument. 

Accuracy: 



I Arguments | Relative Error | 

I I *10**15 I 

I Range | Distribution! RMS | Maximum j 

^^ ^„+ ^ — ^. — ^ + 4 

I ABSCx) < i Uniform 1 0.0638| 0.223 | 
10.25 ^1 III 

I ABSCx) < I Uniform | 0.0913| 0.253 | 
10.95 I III 

i .A — -.--__ X .. ^JL . ^ J 

Error and Exceptional Conditions : 
P I ABS Cx) ^ 1 



ERF^ ERFC (short floating-point r^il) 

Module Name: IHEWEFS 
Entry Points: 

Mathematical PL/ 1 Entry 

Function Name Point 

Error function Cx) ERFCx) IHEEFSF 

Complement of orror ERFCCx) IHEEFSC 

functionCx) 

Function: 

To calculate the error function of x or 
the compl^nent of this function. 



Method 1 



1. 



< X < 1 



Compute ERFCx) by the folloidng 
approximation: 

EHFCx)=x*Ca ♦a-t*x**2+aa*x**«l*. 
**10 



.-»-a-*x 



The coefficients were obtained by the 
minimax approximation in relative error 
of ERFCx)/x as a function of x**2 over 
the range < x**2 < 1. The relative 
error of this approximation is less than 
2**-2«*.6. The value of ERFCCx) is com- 
puted as 

ERFCCx) = 1-ERFCx) 

2. 1 < X < 2.040452 

Compute ERFC Cx) by the following 
approximation : 



ERFCCx) = bo+bj.*2*b2*2**2-i-. 
4- b9*z**9 



where z 



x-Tr 



and To = 1.709472- 



llie 



coefficients were obtained by the minimax 
approximation in absolute error of the 
function f Cz) = ERFCCz+Tq) over the range 
-0.709472 < z < 0.33098. The absolute 
error of this approximation is less than 
2**- 31.5. The limits of this range and 
the value of the origin. Tot were chosen 
to minimize the hexadecimal rounding 
errors. 

Hie value of ERFCCx) within this range is 
between 1/256 and 0.1573. 

The value of ERFCx) is computed as 

ERFCx) = 1-ERFCCx) 

3. 2.040452 < x < 13.306 

compute ERFCCx) by the following 
approximation: 
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ERFCCx) = EXP(~z> ♦F/x 



Where 2 = x**2 and 



F = 



Itie coefficients of F were obtained by 
transforming a roininaax rational approxi- 
mat ion in absolute error of the function 
f(w> = ERFC(x> ♦x«EXPCx**2>over the range 
13.306**-2 < w < 2.040452**-2 Cwhere w=x* 
♦2). This approximation is of the form 



f(w) = 

Cao + aj^*w+a:2*w**2*a3*w**3>/Cbo*bj.*w*w**2) 

The absolute error of this approximation 
is less than 2**-26.1. 

If 2.040452 ^ X < 3.9192, ERF(x) = 1 - 
ERFCCx) 

If 13.306 > K ^ 3.9192, ERFCx) = 1 



4. X > 13.306 

Results 1 and are returned for ERFCx) 
and ERFCCx) respectively. 



5. X < 

Reduce to a case involving a positive 
argument by use of the identities: 

ERFCx) = - ERFC-x) 
and ERFCCx) = 2 - ERFCC-x). 



Effect of an Argument Errors 

The absolute error caused in the result 
is approidroately equal to the absolute 
error in the argianent multiplied by 
EXPC-x**2) . 

ERFCx): As the magnitude of the argument 
increases from 1, the effect of an argu- 
ment error diminishes rapidly- For small 
X, the relative error of the result is of 
the order of the relative error of the 
argument • 

ERFCCx): For x > 1, ERFCCx) is approjci- 
mately EXPC-x**2)/ C2*x) . Thus the rela- 
tive error in the result is approximately 
equal to the relative error in the argu- 
ment multiplied by 2*x**2. For negative, 
or small positive, values of x, the rela- 
tive error in the result is approximately 
equal to the absolute error in the argu- 
ment multiplied by EXPC-x**2). 



Accuracy: 



^ , — - 

1 Arguments 


— ^ . ^ 

1 Relative Error | 


1 




1 *10**6 1 


t— 

1 Range 


1 Distribution! RMS | 


H 

Maximum | 


IHEEFSF 








lABSCx) 


1 Uniform 


1 0.115 I 


0.853 1 


|< 1 

|. 

|1 < 


-+ 


"+ +- 


4 


lABSCx) 


1 Uniform 


1 0.03701 


0.107 1 


|< 2.04 








k 

12.04 < 


-+ " 


-~+ +- 


— — ^-^ 


lABSCx) 


1 Uniform 


1 0.03481 


0.0597 1 


|< 3.9192 








1 ™-„^.™ 


1 


± ± 


J 


IHEEFSC 


1-3.8 < X 


T 

1 Uniform 


T T 

1 0.297 1 


0.941 1 


|< 

k 

|0 < X 


-+- 

1 Uniform 


t 1 

x x 


- - J 


1 0.126 1 


0.692 1 


|< 1 

k 

|1 < X 


-+ 

1 Uniform 


1 t 

~+ +- 

1 0.374 1 


— 4 

1.98 1 


|< 2.04 

J. 

12.04 < X 


-+ 

1 Uniform 


1 1 
1 0.369 f 


4 

1.27 1 


|< 4 

k~ 

|4 < X 


1 Uniform 


1 1 
1 8.22 1 


4 

15.1 1 


|< 13.3 




1 1 




i 


«.x . « 


X X- 


1 



ERF, ERFC Clonq floatinqr- point real) 
Module Name: IHEWEFL 
Entry Points: 

Mathematical PL/I Entry 

Function Name Point 

Error function Cx) ERFCx) IHEEFLF 

Conqploaent of error ERFCCx) IHEEFLC 

fuiKrtionCx) 

Function: 

To calculate the error function of x or 
the complement of this function. 

Method: 

1. S X < 1 

Coaqpute ERFCx) by the follo%fing 
approximation: 

ERFCx) = 

x*Cao-i-aj.«x**2+aa*K**4 + ... .^a, ^♦x**22) 

The <x>efficients were obtained by the 
mini max approximation in relative error 



32 



of ERFCx)/x as a f^onctioe of x**2 over 
the range < x**2 <1. The relatiire 
error of this approximation is less than 
2**-56.9. The value of ERFC is computed 

as 

ERFCCx) = 1-ERFCx) 

2. 1 < X < 2.0404S2 

Compute ERFCCx) by the following 
approximation : 

ERFCCx) = 

bo*bA*z+baZ**2+. . . .♦b 3_e*z**18 

where z = x-^o and To = 1.709472- The 
coefficients were obtained by the minima x 
approxiroation in absolute error of the 
function fCz) = ERFCCz+Tq) over the range 
-0.709472 < z < 0.330948, The absolute 
error of this approximation is less than 
2** -60. 3. The limits of this range and 
the value of the origin, Tq , were chosen 
to minimize the hexadecimal rounding 
errors . 

The value of ERFCCx) within this range is 
between 1/256 and 0.1573. The value of 
ERFCx) is computed as 



3. 



ERFCx) = 1-ERFCCx) 
2.040452 < X < 13.306 



Compute ERFCCx) by the following 
a ppr oximati on : 

ERFCCx) = E)CPC-z)*F/x 

%Aiere z = x**2 and 

F = Co+di/Cz+Cj^) ♦da/Cz+c^)*. •• -^d^/Cz^+c) 

The coefficients of F were obtained by 
transforming a minimax rational approxi- 
mation in absolute error of the function 
f Cw) = ERFCCx)*x*EXPCx**2) over the range 
13.306**-2 < w < 2.040452**-2 Cwhere w=x* 
♦-2) . This approxiroation is of the form 

f Cw) = Cao+ai.*W'i-a2*w**2*. . • .a ♦w**7)/Cb7 
+ bj.*w+b2*w**24'. . . . + b6*w**6+w**7) 

Ihe absolute error of this approximation 
is less than 2**-57.9« 

If 2.040452 < X < 6.092368, IRF(x) = 
1- ERFCCx) 

If 13.360 > X > 6.092368, ERFCx) = 1 

4. X ^ 13.306 

Results 1 and are returned for ERFCx) 
and ERFCCx) respectively. 

5. X < 



Reduce to a case involving positive argu- 
ments by use of the identities: 



ERFCx) = -ERFC~x) 
and ERFCCx) = 2-ERFCC-x) . 

Effect of an Argument Error: 

The ateolute error caused in the result 
is approximately equal to the absolute 
error in the argument iniltiplied by 
EXPC-'X**2). 

ERFCx): hs the magnitude of the argument 
increases from 1, the effect of an argu- 
ment error diminishes rapidly. For small 
X, the relative error of the result is of 
the order of the relative error of the 
argument. 

ERFCCx): For x > 1, ERFCCx) is approxi- 
mately EKPC-x**2)/C2*x). Hi us the rela- 
tive error in the result is approximately 
equal to the relative error in the argu- 
ment multiplied by 2*x**2. For negative, 
or small positive, values of x, the rela- 
tive error in the result is approximately 
equal to the absolute error in the argu- 
ment multiplied by EXPC-x**2). 

Accuracy: 



J.-- . . — . — — _„ 

1 Arguments 


— ^ , ^ .. — -J 

1 Relative Error | 


1 




1 ♦10**15 1 
onf RMS iMaxiirum | 


r— — —' 

1 Range 


iDistributi 


IHEEFLF 






lAESCx)<l 

i 


1 Uniform 
1 


T 0.0257 I 0.193 ] 
1 1 1 


11 < 


lABSCx) 


j Uniform 


1 0.009461 0.02871 


1< 2.04 

i 


1 
1 


1 1 1 

1 1 1 


|2.04 < 


lABSCx) 


1 Uniform 


1 0.008021 0-01391 


|< 6.092 


1 


1 1 ! 


i„ ,„ .,„ .,.,,„. 


1 


± ± _ i 


IHEEFLC 


1-6 < X 


1 Uniform 


1 0.0652 1 0.208 | 


|< 

k— 

|0 < X 


I Uniform 


1 t 1 
1 0.0266 1 0.146 1 


|< 1 
11 < X 


1 Uniform 


1 I 1 
^.+ ^ 4 ^ 

1 0.0913 1 0.426 1 


|< 2.04 
|2.04 < X 


-+ 

1 Uniform 


1 1 1 
1 0.0865 1 0.326 | 


|< 4 




1 1 1 


f 

|4 < X 


1 Uniform 


1 1.96 1 3.51 1 


l< 13.3 




1 1 1 
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FUNCTIONS WITH CXJMPLEX AKGUMENTS 
SQRT (short floating- point complex! 
Module Name: lEEWSQW 

Entry Point: IHESQWO 

Function: 

To calculate the principal value of the 
square root of z, i.e., ~pi/2 < argument 
of result < pi/2. 

Method: 

1. Let SQRTfx+yI> = a+fcl 

2. Let SQRT(CABSCx) + ABS(x+yI))/2) = 
k*SQRT(Wj.*Wa) = C 

v^ = MAX CABS Cx) ,ABS(y)> and 

Va = MIN(ABSCx) ,ABS(y)) 

3. In the special case when either Va = 
or Vj_ » ^2 let 

Wj. = Va and Wa = ^± 

Let k = 1 if v^ = ABSCx) 

k = 1/SQRT(2) if Vi. = ABS(y) 

4. In the general case compute: 

F = SQRT(l/i4-i-Cl/i»)*CVj./Va>**2) 

If ABS(x> is near the underflow thresh- 
old, then take 

WjL = ABSCx), Wa = VjL*2*F, and k = 
1/SQRTC2) 

If v^*F is near the overflow threshold, 
then take 

Wj^ = ABSCx) 74, Wa = Vjl*F/2 and k = 
SQRT C 2) 

In all other cases take 

Wj. = ABSCx) /2, Wa = Vj.*F, pnd k = 1 

5. If c = then a = b = 

If c # and X > 0, then 

a = c, and 
b = y/C2*c) 

if c * and x < 0, then 

a = ABSCy/C2*c)) , and 
b = SIGN(y)*c 

Effect of an Argument Error: 



Let z = r*EXP(hI), and 
SQRTIz) = s*EXPCkI). 

Then the relative error in s is approxi- 
mately half the relative error in r, and 
the relative error in k is approximately 
equal to the relative error in h. 

Accuracy: 

J. , ^ ^-, , ^ ^ 

I Arguments | Relative Error | 
I I ♦lO^^S 1 
^ ^„ ^^ ^ ^4 

I Range | Distribution! RMS | MaxinttUiri | 
I Full range I Exponential | O.SaO ) 2,18 | 

L . JL ^ . X X . J 



SQRT Clonq floating-point complex) 

Module Name: IHEWSQZ 

Entry Point: IHESQZO 

Function: 

To calculate the principal value of the 
square root of z, i.e., -pi/2 < argument 
of result < pi/2. 

Method: 

1. Let SQRTCx+yl) = a+bl 

2. Let SQRTCCABSCx)+ABSCx*yI))/2)= 
k+SQRTCWi+Wa) = c 

Vi = MAXCABSCx) ,ABSCy)) and 
Va = MINCABSCx),ABSCy)) 

3. In the special case when either Va = 
or Vj, » Va let 

Wj. = Va and Wa = v^^ 

Let k = 1 if Vi. = ABSCx) 

k = 1/SQRTC2) if Vj. = ABSCy) 

^. In the general case compute: 

F = SQRTC1/4^C1/1*)*CV:L/Va)**2) 

If ABSCx) is near the underflow thresh- 
old, then take 

Wj. = ABSCx), Wa = Vj,*2*F, and k = 
1/SQRTC2) 

If v^»f« is near the overflow threshold, 
then take 

Wi = ABSCx)/4, Wa = v^^F/2 and k = 
SQRTC2) 
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In all other cases take 

w^ = ftBS{x>/2, W2. = Va*F, and k = 1 

5- If c = then a = b = 

If c ^ and x ^ 0, then 



i^cciiracf 3 



a = c, and 
b = y/C2*c) 

if c # and x < 0, then 

a = ABSCy/C2*c)), and 
b = SIGNCy)*c 



Effect of an Argument Errors 

Let z = r*EXP(hI>, and 
SQRTCz) = s*EXPCkI>. 

Then the relative error in s is approxi- 
mately half the relative error in r, and 
the relative error in k is approximately 
equal to the relative error in h. 



Accuracy: 

J.. . — .^ — -^ . _ — ^«. ^ , ^-^ 

I Arguments | Relative Error | 
I I ♦10**15 I 

J. — ,^^ — ^^^^^ ^„^^„^.^^ ^™„ 4 

I Range (Distribution! RMS | Maximum | 

|. + „^^„^. ^^^^^^ ^^ 

I Full range I Exponential | 0,131 | 0,492 | 

i , ^ X X , -^-L . . ^J 



EXP (short floating-point complex) 

Module Name: IHEWEXW 

Entry Point: IHEEXWO 

Function: To calculate e to the power z- 

Method: 

Let z = X ♦ yl. 



Then REAL ( EXP Cz)) 
and IMAGCEXPfz)) 



EXPCx>*COS(y) 
EXPCx)*SINCy>. 



Effect of an Argument Error: 

Let EXPCx + yl) = s*EXP(kI). 



Then k = y, and the relative error in s 
is approximately equal to the absolute 
error in x. 



Arguments 



,^ _, — , , — ^ 

Relative Error 
♦10*«6 



I Range [Distribution 



|ABS(x) 
|< 170 

|ABS(y> 
l< pi/2 



h 

JABSCx) 
|< 170 

I pi/2 < 
|ABS(y)< 20 



Uniform 



Uniform 



RMS I Maximum 



H 



0.646 



0,628 



2.40 



2.28 



Error and Exceptional Conditions: 

O : ABSCy) > 2**18*pi : error caused in 
real SIN routine CIHEWSNS) 

H : OVERFLOW in real EXP routine 
(IHEWEXS) 



EXP (long floating-point complex) 
Module Name: IHEWEXZ 

Entry Point: IHEEXZO 

Function: To calculate e to the power z. 

Method : 

Let 2 = X + yl . 

Then REALCEXP(z)) = EXP(x)*COSCy) 
and IMAG(EXP(z)) = EXPCx>*SIN(y) . 

Effect of an Argument Error: 

Let EXPCx ♦ yl) = s*EXPCkI). 

Then k = y, and the relative error in s 
is approximately equal to the absolute 

error in x« 

Accuracy: 



Arguments 



--J. — ^ „ . . — ^ 

I Relative Error 
I *10**15 



I 

H- 

I Range (Distribution! RMS | Maximum 

lABSCx) < 1| Uniform | 0,187 | 0.614 

|ABS(y) 1 I I 

l< pi/2 I I I 

lABSCx) I Uniform | 0.200 | 0.819 

|< 20 I I I 

lABSCy) I I I 

l< 20 I I I 
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Error and Exceptional Conditions 3 



O : ABSCyl ^ 2**50*pi : error caused in 
real SIN routine (IHEWSNL) 



H : OVERFLOW in real EXP routine 
flHEWEXL) 



LOG Csbort floating-point complex) 
Module Name: IHEWLNW 

Entry Point: IHELNWO 

Function: 



To calculate the principal value of the 
natural log of 2, i.e., -pi < imaginary 

part of result < pi. 



Method: 

1. Let LOGlx+yl) = a+bl 

2. Then, a = LOG(ABSCx*yl>) and b = 
ATANly^x) 

3. LOG CABS Cx+yl)> is computed as follows: 
Let Vj. = MAXCABS(x>,jyBS(yll and 

Va = l!INCABSCxl,ABSC¥H 

Let t be the exponent of Vj. (i.e., Vj. = 
m*16**t, 1/16 < m < 1) 

Let tj. = t if t < or 

tj. = t-1 if t > and 

s = i6«*tjt 

Then LOG(ABSfx*yI) I = «l*ti*L0GC2> ♦ 
LOGC(Vjl/s)**2 ♦ CVa/s>**2)/2 

Computation of Vj^/s and v^/s are carried 
out by suitable adjustment of the charac- 
teristics of Vj. and Vaj in particular, if 
Va/s « 1, it is taken to be 0. 

Effect of an Argument Error: 

Let z = r*EXPChI) and LOG(z) = u ♦ vl. 

Then the absolute error in u is approxi- 
mately equal to the relative error in r. 
For the absolute error in vC= h = ATANCy, 
x) ) , see corresponding paragraph for 
module IHEWATS. 



Accuracy: 

^ ^ ^ , r«„ ^„«„. ^ 

I Arguments | Relative Error | 
I I ♦lO**^ I 
|« ^ -+ ^^ 4 

I Range f Distribution! RMS | Maximum | 
I Full range I Exponential | 0.396 | 1.89 | 

I ^-A-. J. J. 1 

Error and Exceptional Conditions: 

O: x= y = 0, error in real LOG routine 
(IHEWLNS) 

LOG (long floating- point complex) 

Module Name: IHE^LNZ 

Entry Point: IHELNZO 

Function : 

To calculate the principal value of 
natural log of z, i.e., -pi < imaginary 
part of result < pi. 

Method: 

1. Let LOGCx+yl) = a+bl 

2. Then, a = LOGCABSCx*yII) and b = 
ATAN(y,x) 

3. LOG CABS (x+y I) ) is ccmputed as follows: 
Let Vi = MAX(ABS(x),ABSCy)) and 

Va = MIN(ABSCx> ,ABSCy>) 

Let t be the exponent of Vj. (i.e., Vj. = 
m*16**t, 1/16 < ro < 1> 

Let tjL = t if t < or 

tj. = t-1 if t > and 

s = 16#*ti 

Then LOG CABS Cx+yl>) = 

4*tj.*LOGC2> + LOGCCVj./s)«*2*CVa/s>**2)/2 

Computation of v^/s and Va/s are carried 
out by suitable adjustment of the charac- 
teristics of Vj. and Val in particular, if 
Va/s « 1, it is taken to be 0. 

Effect of an Argument Error: 

Let z = r ♦EXP (hi) and LOGCz) = u ♦ vl. 

Then the absolute error in u is approxi- 
mately equal to the relative error in r. 
For the absolute error in vC= h = ATANCy, 
x)) see the corresponding paragraph for 
module IHEUATL. 
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Accuracy: 



, , ^ — . —J.- ,. — ._ . — -| 

I Arguments | Relative Error | 

I I ♦10**15 I 

^ „^ ^ ^^^^„„„^ .4 

I Range | Distribution! RMS | Maximum f 
1^^ + ^ + 4 

(Full range! Exponential | 0.125 | 0.5i*2 | 

L L «- . X X J 



Error and Exceptional Conditions: 

O : x = y = 0, error caused in log rou- 
tine CIHERLNL) 



The coefficients were obtained by the 
minimax approximation in relative error 
of SINHCx>/x over the range < x**2 < 
0.12011 under the constraint that the 
first term shall be exactly 1. 



u < 



SINHCu) = -SINHC-u). Then 

COSHCu) = SINHCABSCu)) ♦ 1/EXP CABSCu) ) 



Effect of an Argument Error: 



SIN, SINH^ COS, COSH Cshort floating-point 
complex) 



Module Name: 
Entry Points: 



IHEWSNW 



Mathematical 
Function 

Sin z 

Hyperbolic sin 2 
Cos z 
Hyperbolic cos z 

Function : 



PL/ 1 

Name 

SINCz) 

SINHCz) 

COS(z) 

COSH(z> 



Entry 

Point 

IHESNWS 

IHESNWZ 

ih:^nwc 

IHESNWK 



To calculate sin z or hyperbolic sin z, 
or cos z or hyperbolic cos z. 

Method: 

Let z = X + yl. 

Then REALCSIN(z)) = SIN Cx) ♦COSH (y) 
and IMAG(SINCz)) = COS(x) ♦SINHCy) j 

REALfCOSCz)) = COS Cx) ♦COSH Cy) 
and IMAG(COS(z)) = -SIN(x) ♦SINH(y) 1 



REAL(SINH(z)) 
and IMAGCSINHCz)) 



COS(y)^SINH(x) 
SINCy)^COSH{x) I 



REAL(COSHCz)) = COS Cy) ♦COSH Cx) 
and IMAGCCOSHCz)) = SINCy) ♦SINHCx) • 

To avoid making calls to evaluate SINH 
and COSH separately, and thus frequently 
having to evaluate EXP twice for the same 
argiament, SINHCu) is computed as follows: 



1. u > 0.3*165736 
SINHCu) = CEXPCu) 

2. < u < 0.3465736 



l/EXPCu))/2- 



SINHCu)/u is approximated by a polynomial 
of the form a^ ♦ aj,^u+^2 -•• aa*u^^4 C which 
has a relative error of less than 
2^^-26.4). 



Combine the effects on SIN, COS, SINH and 
COSH according to the method of evalua- 
tion described in the above paragraph. 

Accuracy: 



Arguments 



Relative Error 
♦10^^6 



I Range 

L 



I Distribution | 



RMS 



Maximum 



IHESNWS 



lABSCx) 
1^ 10, 
lABSCy) < 



Uniform 



1.17 



3.37 I 



I 

-X- 



IHESNWZ 



lABSCx) I 
|< 10, I 
JABSCy) < 1| 



Uniform 



I 0.878 I 2.75 I 



I 



IHESNWC 

r — -" 

JABSCx) 
|< 10, 
JABSCy) 



I Uniform [ 1.18 | 3.23 | 



I 



I 

-X- 



I 






IHESNWK 



JABSCx) 
|< 10, 
JABSCy) 
|< 1 



I Uniform | 0.968 | 3.11 | 
I III 



I 

-X- 



I 

-X- 



Error and Exceptional Conditions: 

O : IHESNWS, IHESNWC: 

ABSCx) ^ 2^^18+pi : error caused in 
real SIN routine CIHEWSNS) 

IHESNWZ , IHESNWK : 

ABSCy) ^ 2^^18+pi : error caused in 
real SIN routine CIHEWSNS) 

H J OVERFLOW in real EXP routine 
CIHEWEXS) 
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SIN, SINK. COS. COSH (long floating-point 
complex) 



Module Name: 



Entry Points: 



IHEIiSNZ 



Mathematical 




PL/I 


Entry 


Function 




Name 


Point 


Sin z 




SINCz) 


IHESNZS 


Hyperbolic sin 


z 


SINH(z) 


IHESNZZ 


Cos z 




COSCz) 


IHESNZC 


Hyperbolic cos 


z 


COSH(z) 


IHESNZK 



Function: 

To calculate sin z or hyperbolic sin z^ 
or cos z or hyperbolic cos z. 



Method: 

Let z = X + yl. 

Then REAL(SINCz)) 
and IMAG(SIN(z)) 



and 



REALCCOS(z)) 
IMAGCCOS(z)) 



= SIN (x> ♦COSH Cy> 
= COS(x>*SINH(y)| 

= COSCx)*COSHCy) 

= -SINCx>*SINH(y)| 



REAL(SINHCz>) = COS Cy) ♦SINHCx) 
and IMAGCSINHCz)) = SINCy) *COSH(x) | 

REiUlitCOSHCz)) = COS(y)*COSH(x> 
and IMAG(COSH(z)> = SIN(y> ♦SINHCx) • 

To avoid making calls to evaluate SINH 
and COSH separately, and thus frequently 
having to evaluate EXP twice for the same 
argument, SINHCu) is computed as follows: 



1. u > 0.481212 

SINHCu) = CEXPCu) 

2. < u < 0.481212 



l/EXP(u))/2 



Accuracy s 

J. «. , , ^ ^ ^ 

I Arguments | Relative Error | 

I I ♦lO+^lS I 

^ ^„^ + ^ 4 

I Range | Distribution! RMS j Maximum | 

I ^ ^ i. . — i , ^j. J 

IHESNZS 

^ ^ , ^«, ^ , ^ 

JABSCx) I Uniform | 2.01 | 113 | 

N 10» I II I 

|ABS(y) < 1| I I I 

IHESNZZ 

lABSCx) I Uniform | 0.229 | 0.6411 

1^ lOi I I I I 

|ABS(y) < 1| II I 

IHESNZC 

J. . — -^ ^ ^ ^ ^ 

lABSCx) I Uniform | 0.311 | 3.83 | 

1^ 10, I III 

|ABS(y) < 1| II I 

L ^ ^J. J. JL J 

IHESNZK 

|ABS(x) I Uniform | 0.250 | 0.730 | 

l^ 10, I II I 

lABSfy) I III 

1^ 1 I II I 

L ^ 1. , JL .. X J 



Error and Exceptional Conditions: 

O : IHESNZS, IHESNZC: 

ABSCx) ^ 2**50*pi: error caused in 
real SIN routine (IHEWSNL) 

IHESNZZ, IHESNZK: 

ABSCy) ^ 2**50*pi: error caused in 
real SIN routine C IHEWSNL) 

H : OVERFLOW in real EXP routine 
CIHEWEXL) 



SINH(u)/u is approximated by a polynomial 
of the fifth degree in n**2 which has a 
relative error of less than 2**-56.07 

The coefficients were obtained by the 
minima X approximation in relative error 
of SINHCx)/x over the range < x**2 < 
0.23156 under the constraint that the 
first term shall be exactly 1. 

3. u < 

SINHCu) = -SINHC-ul. Then 

COSHCu) = SINH(ABSCu)) ♦ l/EXPCABSCu) I . 

Effect of an Argument Error: 

Combine the effects on SIN, COS, SINH and 
COSH according to the method of evalua- 
tion described in the above paragraph. 



TAN, TANH ishort floating-point complex) 
Module Name: IHEWTNW 

Entry Points: 

Mathematical PL/I Entry 

Function Name Point 

Tan z TANCz) IHETNWN 

Hyperbolic tan z TANH(z) IHETNWH 

Ftinction : 

To calculate tan z or hyperbolic tan z. 
Method : 

Let z = X ♦ yl. 
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Then REJULCTiyilCzl) = 

TAN(x>*Cl - TANHCy)**2}/ 

CI ♦ CTMi(x)*TANH(y))**2>t 



IMAGCTANCz)) = 

TANHCy)*(l + TANCx)**2)/ 

CI ■*• (TAN(xl*TANeCy)l**2). 

TANHCz) = - CTAN{2l))I, 



Effect of an Argument Error: 

The absolute error caused in the result 
is approxiroately equal to the absolute 
error in the argument divided by ABSCCOS 
Cz)**2) for IHETNWN, or divided by ABS 
CCOSHCz>**2> for IHETNWH. The relative 
error caused in the result is approxi- 
mately twice the absolute error in the 
argument divided by ABS (SIN C2*z) > for 
IHETNWN, or divided by ABS(SINHC2*z) ) for 
IHETNWH. 



Accuracy: 

I 

Arguments 



. ^ , — ^^ — . ^--^ 

I Relative Error | 
I I ♦10**6 I 

^^^^ — ^ — ^ „ — ^.-.„„„^ ^ 

I Range | Distribution! RMS | l^ximum | 

L ^^ JL- ._ J. ^-1.. J 

IHETNWN 



|ABS(x) < 1| Uniform 
|ABSCy> < 9| 

L . X - ._ 



I 0.532 I 
I I 



2-87 I 



IHETNWH 



J.-, ^- 

jABS(x) < 9J 
|ABS(y) < 1| 

L ^J.- 



Uniform 



I 0,524 I 



2.71* I 

I 

J 



Error and Exceptional Conditions: 

I : OVERFLOW 

O : ABSCu) > 2**18*pi, where 
u = X for IHETNWN, 
u = y for IHETNWH, 

H : OVERFLOW or ZERODIVIDE in real TAN 
routine CIHEWTNS) 

TAN^ TANH Clonq floating-point complex) 
Module Name: IHEWTN2 
Entry Points : 



Mathematical 
Function 

Tan 2 
Hyperbolic tan z 



PL/I 
Name 
TANCz) 
TANHCz) 



Entry 
Point 

IHETNZN 
IHETNZH 



Function: 

To calculate tan z or hyperbolic tan z- 

Method : 

Let z = X ♦ yl. 

Then REALCTANCzH = 

TANCx)*Cl - TANH(y)**2>/ 
CI ♦ CTAN(x>*TANHCy>)**2>, 

IMAGCTANCzH = 

TANH(y)*Cl ♦ TANCx)**2)/ 
(1 + (TANCx)*TANHCy))**2). 

TANe(z> = - CTAN(zI>)I. 

Accuracy: 



Arguments 



I Relative Error 
I ♦lO^+lS 



I Range 



i Distribution! RMS j Maximum 



IHETNZN 

lABSCx) < 1| Uniform 
JABSCy) < 9 1 



I 0.172 I 0.709 I 



I 

-J 



I 



IHETNZH 



0.692 



J.-. -^ -^ .^- 

|ABS{x) < 9 1 Uniform | 0.174 | 
lABS(y) < 1| I I 

I - . X - X X- 

Error and EXCEPTIONAL Conditions: 

I : OVERFLOW 

O : ABSCu) > 2**50*pi, where 
u = X for IHETNZN, 
u = y for IHETNZH. 



OVERFLOW or ZERODIVIDE in real TAN 
routine CIHEWTNL) 



JITANi, ATANH Cshort floating-point complex) 

Module Name: IHEWATW 

Entry Points: 

Mathematical PL/I Entry 

Function Name Point 

Arctan z ATANCz) IHEATWN 

Hyperbolic arctan z ATANH (z) IHEATWH 

Function: 

To calculate arctan z or hyperbolic arc- 
tan z. 

Method: 

Let z 



X + yl. 
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Then REALlATMIH(z)) = 



IMAGIATANHCz)) = 



CATANHC2*3C/ 

(ATMlC2*y, 
Cl-x*x-y*y)))/2 



and ATANfz) 



= -CATANHCzI)}!, 



Effect of an Argument: Error: 

The absolate error in the result is ap- 
proximately equal to the absolute error 
in the argument divided by CI ♦ z**2) in 
the case of IHEATWN, or by (1 - z**2) in 
the case of IHEATWH, Thus the effect may 
be considerable in the vicinity of z = 
± II CIHEATWN) or ± 1 (IHEATWHK 

Accuracy: 



^- ^ , ^ ^_^ , ^ 

I Arguments | Relative Error | 
I I ♦lO^^e I 

I Range | Distribution) RBS 1 Maximum j 

L^, , . J. , X JL , I 

IHEATWN 



1.05 I 

^-J 



I I II 

|Full range I Exponential | 0.205 | 

I . X X , X- 

IHEATWH 

J. ^ ^ ^ ^ 

I Full range I Exponential | 0.224 | 1. 22 ( 

L : X , X X-. J 



Error and Exceptional Conditions: 



P : IHEATWN: 
IHEATWH: 



z = ±11 
z = ±1 



ATAN^ ATANH (long floating-point ccwnplex) 



Module Name: 
Entry Points ; 



IHEWATZ 



Mathematical 
Function 
Arctan z 
Hyperbolic arctan z 



PL/ 1 
Name 
ATAN(z) 
ATANH(z) 



Entry 

Point 

IHEATZN 

IHEATZH 



Function: 



To calculate arctan z or hyperbolic arc- 
tan z. 



Method : 

Let z = X -1^ yl. 

Then REAL(ATANHCz) > = CATANH(2*x/ 

(H^x*x+y*y>>)/2 

(ATAN(2*y, 
(l-x*x-y*y)))/2 



IMAGCATANHCz)} = 



and ATAN(z) 



= -CATANH(zI))I. 



Effect of an Argument Error: 

The absolute error in the result is ap- 
proximately equal to the absolute error 
in the argument divided by (1 ♦ z**2) in 
the case of IHEATZN, or by (1 - z**2) in 
the case of IHEATZH. Thus the effect may 
be considerable in the vicinity of z = 
± I (IHEATZN) or ± 1 (IHEATZH). 



Accuracy: 

I ^, ^ _ ^ 

I Arguments [ Relative Error | 
I I *10**15 I 

^ r ~i T—^ i 

I Range | Distribution! RMS | Maximum | 

I ^ X X X J 

IHEATZN 

J.— — ^ J. ^. — ^ ^^ ^ ^ 

[Full range I Exponential | 0.0517[ 0.^*38 | 

I . X , X -X J 

IHEATZH 

I . ^ ,- ^ -^ 

I Full range (Exponential | 0.0562| 0.4*09 | 

i. . X- X„, X . J 

Error and Exceptional Conditions: 

P s IHEATZN: z = ±11 

IHEATZH: z = ±1 



^♦0 



CHAPTER ^: ABRAY INDEXERS AND BUILT-IN FUNCTIONS 



The Library supports the array built-in 
functions SUM, PROD, POL¥, ALL and ANY, and 
also provides indexing routines for han- 
dling simple (i.e., consecutively stored) 
and interleaved arrays. 

Input Data 

The array function modules are distin- 
guished from the other Library modules in 
that they all accept array arguments ana 
perform their own indexing, whereas the 
other modules require that indexing should 
De handled by compiled code. Calls to con- 
version routines are included in the SUM, 
PROD and POLY modules with fixed- point 
arguments, so that these arguments are con- 
verted to floating-point as they are 
accessed (it should be noted that it is a 
requirement of the language that the 
results from these modules be in floating- 
point) . On the other hand, the conversions 
necessary for the ALL and ANY modules (the 
arguments must be converted to bit string 
arrays) are not part of the modules and 
must be carried out before the modules are 
invoked. 

Any restrictions on the admissibility of 
arguments are noted under the headings 
'Range' and 'Error and Exceptional 
Conditions" . 

Range : This states any ranges of arguments 
for which a module is valid. Arguments 
outside the ranges given are assumed to 
have been excluded before the module is 
called. 

Error and Exceptional Conditions ; These 
cover conditions which may result from the 
use of a routine! they are listed in four 

categories: 

P — Programmed conditions in the module 
concerned. Programmed tests are 



made where this is not too costly 
and, if an in'^^lid argument is 
found, a branch is taken to the 
entry point IHEERRC of the execution 
error package C EXEP> . This results 
in the printing of an appropriate 
message and in the ERROR condition 
being raised. 



I — Interruption conditions in the 

module concerned. For those rou- 
tines where SIZE and FIXEDOVERFLOW 
are detected by programmed tests or 
where hardware interruptions may 
occur, the OVERFLOW, UNDERFLOW, and 
(when the conversion package is 
called) SIZE conditions pass to the 
ON condition error handler (IHEWERR) 
and are treated in the normal way. 
The machine is assumed to be enabled 
for all interruptions except signi- 
ficance, which is masked off. 



O - — Programmed conditions in modules 
called by the module concerned. 
These occur when invalid arguments 
are detected in the module called. 

H ~ As I, but the interrupt conditions 
occur in the modules called by the 
module concerned - 



Effect of Hexadecimal Truncation 

See the corresponding section in the 
introduction to Chapter 3 for guidance to 
the accuracy of SUM, PROD, and POLY. If 
fixed- point arguments are passed to these 
functions, further errors may be introduced 
by conversions. 

A summary of the Library array modules 
is given in Figures 1 and 8. 



J. — . .^ .^-.— . — ». . ^, , — , , ^ — - — . -.. .-.-. . ^ ^ 

I I Simple arrays, and | Interleaved string | 
I I interleaved arrays of | arrays with fixed- | 
I I variable- length strings | length elements | 
^„„ ^+ ^ — ^ — , . ^„.+^.„ „. — ^^^^ ^-„»^4 

I Indexers | IHEWJXS | IHEWJXI | 

I ALL, ANY I IHEWNLl j IHEWNL2 | 

|.^^ „^j. — .^^ ^„^^ ^ j.„ ^„„ ^ 4 

t Note : IHEWJXI is used for indexing through interleaved arithmetic arrays. [ 

L . __^«, . .«_ . « . «_ ^ . J 

Figure 7. Bit String Array Functions and Array Indexers 
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PL/I 
function 



SUM real 

complex 



. — ^ ^-. — ^ . , ^. ^ 

Fixed -point | Floating-point argiiinents | 

arguments |- -i— ^ ^ -| 

I Short precision | Long precision | 

Simple I Interleaved! Simple | Interleaved} Simple (Interleaved | 

IHEWSSFI IHEWSMF ( IHEWSSG | IHEWSMG | IH£WSSH| IHEWSMH | 
IHEWSSXJ IHEWSMX | IHEWSSG J IHEWSWG | IHEWSSH | IHEWSMH \ 



PROD real 

complex 

POLY real 

complex 



f 



Figure 8 , 



I IHEWPSSJ IHEWPDS | IHEWPSL| IHEWPDL | 
I IHEWPSWI IHEWPDW | IHEWPSZt IHEWPDZ | 
+ ± ^+ J ^ 

I IHEWYGS t IHEWYGL j 
I IHEWYGW I IHEWYGZ | 
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IHEWPSFf IHEWPDF 
IHEWPSXl IHEWPDX 

X 

IHEWYGF 
IHEWYGX 



ABJRAY INDEXERS 
Index er for Simple Arrays 
Module Name: IHEWJXS 
Entry Points: 



Element 

Address 
Bit addresses 
Byte addresses 



Entry 

Point 

IHEJXSI 

IHEJXSY 



Function: 

To find the first and last elements of an 
array. "Iheir addresses are returned, in 
general registers and 1 respectively, 
as bit addresses (IHEJXSI) or byte 
addresses (IHEJXSY) . 

Method: 

The address of the virtual origin B of 
the array (i.e., the address that would 
correspond to the element A(0,..0)) is 
obtained as a byte address for IHEJXSY, 
or a bit address for IHEJXSI, by refer- 
ring to the first word of the array dope 
vector (ADV). 



Address of first element = B ♦^ ^±^i 

i=l 
n 

Address of last element = ^ ♦ J^ ^i^i 

i=l 

where M is the mtiltiplier for the ith 
dimension 

L is the lower bound for the ith 
dimension 

U is the upper bound for the ith 
dimension, and 

n is the number of dimensions. 



Range: 

< nuinber of dimensions < 2 ♦♦22 
Indexer for Interleaved Arrays 
Module Name: IHEWJXI 
Entry Points: 









Entry 


Operation 






Point 


In it ial i za t ion 


for 


bit 


IHEJXII 


addresses 








In it ial iza tion 


for 


byte 


IHEJXIY 


addresses 








Elements after 


the 


first 


IHEJXIA 



Function: 

To find the next element of an array and 
to return its bit or byte address in gen- 
eral register 1. 

Entry point IHEJXII is used to initialize 
the routine for bit addresses and to pro- 
vide the address of the first elem^it in 
the array? IHEJXIY does the same for byte 
addresses. Entry point IHEJXIA is used 
thereafter to obtain the addresses of 
subsequent elements of the arrayi one 
address is returned for each entry into 
the routine. 

Method: 

Arrays are stored in row major order. 
Let Li be the lower bound and O^ the 
upper bound of the ith dimension, and n 
the number of dimensions. Starting with 
the element A(Lj.,La, - . • • • -L^) , the rou- 
tine varies the subscripts through their 
ranges to A(UjL,Ua,. . . • ..Un) , changing the 
nth subscript mc^t rapidly? in this way 
the elements are referenced in the order 
in which they are stored. 

*nie routine does not deal with actual 
subscript values but calculates the 
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extient. E£C= 



^i ■" ^ i * 1^ ^f each dimen- 



sion and uses this as a cotint that varies 
from E|_ to 1 for subscript values Li to 
B I, h •base address" for ea€±t dimension 
is maintained and, for the ith dimension, 
is defined as the address of the elaoaent 
with ith subscript equal to its lowest 
bound Li and with all other subscripts at 
their current values. 

Thus initially the base addresses are all 
equal to the address of ACLi.,L2# • • • • • • 

Ln) • Each subsequent element address is 
generated from the previous one by adding 
the multiplier M,^ from the array dope 
vector iM)V} and reducing the subscript 
count by 1« When the count for the ith 
dimension has been reduced from E^ to 1 
it is reset to Ej_, ^±--1 is added to the 
Ci-l)th dimension's base address and the 
count for this dimension is decreased by 
one. 

This new base is the starting point for 
further increments by Mn» When a new 

base address is calculated, the base 
addresses for all higher dimensions 

C(i+1), Ci-i'2) ....... n) is set equal to 

the ith base address. 

Range : 

< number of dimensions < 2**22 



ARRAY itJNCTIONS 

ALL CX), ANY CX) 
Module Names : 



modules) • The result is a scalar bit 
string of length equal to the greatest of 
the current lengths of the elements of X. 



ALLCX): the ith bit of the result is 1 if 
the ith bits of all the elements of X 
exist and are 1; otheriid.se it is 0. 

A'NYCXI: the ith bit of the result is 1 if 
the ith bit of any element of X exists 

and is 1| otherwise it is 0. 

Method: 

For byte-aligned string arrays, AND 
CIHEMBSA} and OR (IHEWBSO) are used for 
ALL and ANY respectively! for string 
arrays with any alignment BOOL (IHSIBSF) 
is used with appropriate parameter bits. 

The eleoients of the array are passed to 
IHEWBSA, IHEWBSO, or IHEWBSF one at a 
time, and the result is developed in the 
target field- For the first call to any 
of these logical modules the first ele- 
ment of the array serves as both first 
and second source arguments. For subse- 
quent calls, the result already developed 
in the target field is the first argument 
and the next element of the array is the 
second argument. 

Range: 

Bit strings are limited to a maxim\im of 

32,767 bits. 



SUM CX) 



Arguments 
Simple arrays and interleaved 
arrays with variable- length 
elements 

Interleaved arrays with fixed- 
length elements 

Entry Points: 



PL/I Function 
ALLCX), ANYCX), byte-aligned 

ALL(X), any alignment 

ANYCX), any alignment 

Function: 



Module 
Name 
IHEWILl 



IHEWNL2 



Entry 

Point 

IHENLIA 

IHENL2A 

IHENLIL 
IHENL2L 

IHENLIN 
IHELN2N 



Module Names and Entry Points: 

Simple Arrays 

Modiile Entry 

Arguments Name Point 

Fixed, real IHEWSSF IHESSFO 

Fixed, complex IHEWSSX IHESSXO 

Short float 

real IHEWSSG IHESSGR 

complex IHEWSSG IHESSGC 

Long float 

real IHEWSSH IHESSHR 

complex IHEWSSH IHESSHC 

Interleaved Arrays 



Hie argument X is a bit string array Cany 
necessary conversion having been per- 
formed prior to the invocation of these 



Arguments 

Fixed, real 

Fixed, complex 

Short float 

real 
complex 

Long float 

real 
complex 



Module 

Name 
IHEWSMF 
IHlSfSMX 



IHEWSMG 
IHEWSMG 

IHEWSMH 
IHEWSMH 



Entry 

Point 

IHESMFO 

IHESMXO 

IHESMGR 
IHESMGC 

IHESMHR 
IHESMHC 
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Function: 



Method : 



To produce a scalar with a value which is 
the sum of all the elements of the array 
argument. 

Method: 

The elements of the array are added to 
the curra:it sum in row major order. 

For fixed -point arguments each element is 
converted to floating-point by using the 
PL/I Library conversion package. 

For a complex argument, the summation of 
the real parts is performed before the 
summation of the imaginary parts is begun 
in modules IHEWSSG and lEEWSSH, while the 
two sums are developed concurrently in 
other modules. 

Error and Exceptional Conditions : 

I : OVERFLOW, ONDERFL0W 

e : IHEWSSF, IHEWSSX, IHEWSMF, IHEWSMX: 
ABS (element of the arrayl > 7.2*10** 
75: SIZE condition caused in conver- 
sion package 



PROD (X) 



Module Names and Entry Points: 
Simple Arrays 





Module 


Entry 


Arguments 


Name 


Point 


Fixed, real 


IHBIPSF 


IHEPSFO 


Fixed complex 


IHEWPSX 


IHEPSXO 


Short float 






real 


IHEWPSS 


IHEPSSO 


complex 


IHEWPSW 


IHEPSWO 


Long float 






real 


IBBIPSL 


IHKPSLO 


complex 


IHEWPSZ 


IHEPSZO 


Interleaved Arrays 




Module 


Entry 


Arguments 


Name 


Point 


Fixed, real 


IHEWPDF 


IHEPDFO 


Fixed, complex 


IBIWPDX 


IHEPDXO 


aiort float 






real 


IHjBIP DS 


IHEPDSO 


complex 


IHEWPDW 


IHEPDWO 


Long float 






real 


IHEWPDL 


IHEPDLO 


complex 


IHB»PDZ 


IBEPDZO 


Function: 







ThB elements of the array are used in row 
major order to multiply the current 
product. 

For fixed- point arguments, each element 
is converted to floating-point by using 
the PL/I Library conversion package- 

Error and Exceptional Conditions: 

I : OVERFLOW, UNDERFLOW 

H : IHEWPSF, IHEWPSX, IHEWPDF, IHEWPDX: 
ABSCelement of the array) > 7.2*10** 
75: SIZE condition caused in conver- 
sion package 

POLY (A,X) 

Module Names and Entry Points: 







Module 


Entry 


Arguments 




Name 


Point 


Fixed, real 








vector X 




IHEWYGF 


IHEYGFV 


scalar X 




IHEWYGF 


IHEYGFS 


Fixed, complex 






vector X 




IHEWYGX 


IBEYGXV 


scalar X 




IHEWYGX 


IHEYGXS 


Short float. 


real 






vector X 




IHEWYGS 


IHEYC^V 


scalar X 




IHEWYGS 


IHEYGSS 


Short float. 


complex 




vector X 




IHEWYGW 


IBEYGWV 


scalar X 




IHEWYGW 


IHEYGWS 


Long float, real 






vector X 




IHEWYGL 


IHEYGLV 


scalar X 




IHEWYGL 


IHEYGLS 


Long float, compl< 


ex 




vector X 




IHEWYGZ 


IHEYGZV 


scalar X 




IHEWYGZ 


IHEYGZS 


mction: 








Vector X: 









To produce a scalar with a value %rtiich is 
the product of all the elements in the 
array argument. 



Let the arguments be arrays declared as 
A(msn) and X(p:g). Then the function 
computed is: 

n^ro j-1 

A(m) ♦ J2 A(m ♦ j)* n X (p ♦ i) 

j=l i=0 

unless n = m, when result is A Cm) . 

Ifg-p<n-m-l, then, for p ♦ 
i > g, X(p ♦ i) = X(g). 



Scalar X: 

This may be interpreted as a special case 
of vector X, that is, a vector with one 
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elament, XC1>, which is equa.l to X. Then 
the function computed is: 

J] him ♦ j)*X**j 

h floating-point result is obtained in 
both cases. 

Method: 

1. Vector X, (q-p>n--i!i-l): 

POLY(A,X> is evaluated by nested multi- 
plication and addition, i.e., 

C...CACn)*XCk) + ACn-l))*XCk-l) + 

ACn-2))* .•. ^' ACnH-1) ) *X(p> + him} 

where k=p"i-n~ni- 1. 

2. Vector X, {q - p < n - m - 1) : 

In the expression above, the terms in X 
mth subscript ranging from k down to q 



are all raade «3ual to X(q). The evalua- 
tion is treated as for scalar X until 
sufficient terms in X have been iiede 
equal to XCq), when the computation con- 
tinues as in Cl.>. 



3. Scalar X: 

Terms in X with subscript longing from k 
to p are equal to X. 

For fixed-point arguments each element is 
converted to floating-point, by using the 
PL/I Library conversion package. 

Error and Exceptional Conditions: 

I : OVERFI/)W, UNDERFLOW 

H : IHEW¥GF, IHEWYGX: 

ABS( element of the array) > 7.2*10** 
75: SIZE condition caused in conver- 
sion package 
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INDEX 



g •and* operator 3 
I I cx>n catenate operator 
-• 'not* operator 3 
I 'or* operator 3 



ABS 

complex fixed- point 15 

complex floating- point 15 
absolute error, definition 19 
accuracy, in 

aritliroetic operations 8 

raatheroatical functions 17 
ADD 

complex arguments 13 

real arguments 13 
algorithms 18 
ALL array function 45 
'and' operator 3 
ANY array function 45 
array functions 45 
array indexers 

interleaved arrays 44 

simple arrays 44 
assignment operations 

bit string 4,5 

character string 7 
AT AN 

complex arguments 41 

real arguments 27 
ATAND (real arguments) 27 
ATANH 

complex arguments 41 

real arguments 31,32 



bit string operations 3-6 

and 3 

assign, general 4 

assign/fill 5 

BOOL 6 

comparison, byte aligned 

comparison, general 4 

concatenate 4 

INDEX 5 

not 3 

or 3 

REPEAT 4 

SUBSTR 5 
BOOL 6 
built-in functions 

arithmetic 8 

array 42 

bit string 5 

character string 7 

mathematical 17 



character string operations 
assign/fill 7 
compare 6 
con ca t Qia t e 6 



INDEX 7 

REPEAT 6 

SOBSTR 7 
comparison operator 

bit strings 4 

character strings 6 
concatenation operator 

bit string 4 

character string 6 
QOS 

complex arguments 38,39 

real arguments 23,24 
COSD (real arguments) 23,24 
COSH 

complex arguments 38,39 

real arguments 28,29 

DIVIDE (complex fixed- point) 14 

division operator 

complex fixed-point 10,11 
complex floating point 11 



ERF (real arguments) 32,33 
ERFC (real arguments) 32,33 
error conditions, in 

arithmetic operations and functions 
array indexers and functions 17 
mathematical functions 43 
exponentiation operator 
ccanplex operations 

floating-point deponents 12 
integer exponents 12 
real operations 

floating-point exponents 9 
integer exponents 9 



fill operations 
bit string 5 
character string 



HIGH 16-17 



IHEABU 

see ABS (complex fixed-point) 
IHEABV 

see ABS (complex fixed- point) 
IHEABH 

see ABS (complex floating-point) 
IHEABZ 

see ABS (complex floating-point) 
IHEADD 

see ADD (real arguments) 
IHEADV 

see ADD (complex arguments) 
IHEAPD 

see shift-and-assign, shift -and- load 
(real operations) 
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IHEATL 

see ATAN (real arguroentsl ; ATAND (real 
arguments) 
IHEATS 

see ATAN (real arguments)! ATAND (real 
arguit^nts) 
IHEATW 

see ATAN C complex arguments); ATANH 
(complex arguments) 
IHEATZ 

see AT AN (complex arguments) j ATANH 
(complex arguments) 
IHEBSA 

see "and* operator 
IHEBSC 

see comparison operator (bit string, 
byte -aligned) 
IHEBSD 

see comparison operator (bit string, 
general) 
IHEBSF 

see BOOL 
IHEBSI 

see INDIX (bit string) 
IHEBSK 

see concatenation operator (bit string) i 
REPEAT (bit string) 
IHEBSM 

see assignment operations (bit string) j 
fill operations (bit string) 
IHEBSN 

see 'not" operator 
IHEBSO 

see 'or* operator 
IHEBSS 

see SUBSTR (bit string) 
IlECSC 

see comparison operator (character 
string) 
IHECSI 

see INDEX (character string) 
IHECSK 

see concatenation operator (character 
string) ; REPEAT (character string) 
IHECSM 

see assignment operations (character 
string) ! fill operatioi^ (character 
string); HIGH; LOW 
IHECSS 

see SUBSTR (character string) 
IHEDVU 

see DIVIDE (complex fixed-point) 
IHEDW 

see DIVIDE (complex fixed-point) 
IHEDZW 

see division operator (complex 
floating-point) 
IHEDZZ 

see division operator (complex 
floating-point) 
IHEEFL 

see ERF (real argxjments) ; ERFC (real 
arguments) 
IHEEFS 

see ERF (real arguments); ERFC (real 
arguments) 
IHEEXL 

see EXP (real arguments) 



IHEEXS 

see EXP (real arguments) 
IHEEXW 

see EXP (complex arguments) 
IHEEXZ 

see EXP (complex arguments) 
IHEKTL 

see ATANH (real arguments) 
IHEHTS 

see ATANH (real arguments) 
IHEJXI 

see array indexers (interleaved arrays) 
IBEJXS 

see array indexers (simple arrays) 
IHELNL 

see LOG (real arguments); LOG2| LOGIO 
IHELNS 

see LOG (real arguments); L0G2; LOGIO 
IHELNW 

see LOG (complex arguments) 
IHELNZ 

see LOG (complex arguments) 
IHEMPU 

see MULTIPLY (complex fixed-point) 
IHEMPV 

see MULTIPLY (complex fixed-point) 
IHEMXB 

see MAX (real arguments); MIN (real 
arguments) 
IHEMXD 

see MAX (real arguments); MIN (real 
arguments) 
IHEMXL 

see MAX (real arguments); MIN (real 
arguments) 
IHEMXS 

see MAX (real arguments); MIN (real 
arguments) 
IHEMZU 

see multiplication operator (complex 
fixed- point) ; division operator (com- 
plex fixed-point) 
IHEMZV 

see multiplication operator (complex 
fixed-point); division operator (com- 
plex fixed- point) 
IHEMZW 

see multiplication operator (complex 
floating-point ) 
IHEMZZ 

see multiplication operator (complex 
f loati ng- poi nt ) 
IHENLl 

see ALL; AblY 
IHENL2 

see ALL, M^ 
IHEPDF 

see PROD 
IHEPDL 

see PROD 
IHEPDS 

see PROD 
IHEPDW 

see PROD 
IHEPDX 

see PROD 
IHEPDZ 

see PROD 
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I HEPS F 

see PROD 
IHEPSL 

see PROD 
IHHPSS 

see PROD 
IHEPSW 

see PROD 
IHEPSX 

see PROD 
IHEPSZ 

see PROD 
IHESHL 

see SINH (real argiaments); OOSH (real 
arguments) 
IHESHS 

see SINH (real arguments); COSH (real 
arguments ) 
IHESMF 

see SUM 
IHESMG 

see SOM 
IHESMH 

see SUM 
IHESMX 

see SOW 
IHESNL 

see SIN (real arguments); SIND (real 
arguments) I COS (real arguments); COSD 
(real arguments) 
IHEa^S 

see SIN (real arguments); SIND (real 
arguments); COS (real arguments); COSD 
(real arguments) 
IHESNW 

see SIN (complex arguments); SINH (com- 
plex arguments) ; COS (complex argu- 
ments); COSH (complex arguments) 
IHESNZ 

see SIN (complex arguments); SINH (com- 
plex arguments) i C<3S (complex argu- 
ments); COSH (complex arguments) 
IHESQL 

see SQRT (real arguicents) 
IHESQS 

see SQRT (real arguments) 
IHESQW 

see SQRT (ccxnplex arguments) 
IHESQZ 

see SQRT (complex arguments) 
IHESSF 

see SUM 
IHESSG 

see SUM 
IHESSH 

see SUM 
IHESSX 

see SUM 
IHETHL 

see TANH (real arguments) 
IHETHS 

see TANH (real arguments) 
IHETNL 

see TAN (real arguments); TAND (real 
arguments ) 
IHETNS 

see TAN (real arguments); TAND (real 
arguments) 



IHETNW 

see TAN (complex arguments); TANH (com- 
plex argiwients) 
IHETNZ 

see TAN (complex arguments); TANH (com- 
plex arguments) 
IHEXIB 

see exponentiation operator (real opera- 
tions, integer exponents) 
IHEXID 

see exponentiation operator (real opera- 
tions, integer exponents) 
IHEXIL 

see exponentiation operator (real opera- 
tioi^, integer exponents) 
IHEXIS 

see exponentiation operator (real opera- 
tions, integer exponents) 
IHEXIU 

see exponentiation operator (complex 
operations, integer expon^its) 
IHEXIV 

see exponentiation operator (complex 
operations, integer exponents) 
IHEXIW 

see exponentiation operator (complex 
operations, integer exponents) 
IHEXIZ 

see exponentiation operator (complex 
operations, integer exponents) 
IHEXXL 

see exponentiation operator (real opera- 
tions, floating-point exponents) 
IHEXXS 

see exponentiation operator (real opera- 
tions, floating-point exponents) 
IHEXXW 

see exponentiation operator (complex 
operatiox^, floating-point exponents) 
IHEXXZ 

see exponentiation operator (complex 
operations, floating-point exponents) 
IHEYGF 

see POLY 
IHEYGL 

see POLY 
IHEYGS 

see POLY 
IHEYGW 

see POLY 
IHEYGX 

see POLY 
IHEYGZ 

see POLY 
INDEX 

bit string 5 
character string 7 
indexers 

see array indexers 



LOG 

complex arguments 37 
real argmnents 21,22 

LOG2 (real arguments) 21,22 

LOGIO (real argixnents) 21,22 

LOW 7 
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mathonatical functions 19 
MAX Creal arguments) 13 
MIN (real arguments) 13 
Modules Names 1, 2 
multiplication operator 

complex fixed-point 10,11 
complex floating-point 11 
MULTIPLY (complex fixed-point) 



l«l 



shift- and- assign^ shift-and-load (real 

operations) 10 
SIN 

complex arguments 38,39 

real arguments 23, 2i* 
SIND (real arguments) 23,2*1 
SINH 

complex arguments 38,39 

real arguments 38,39 



• not • operator 3 



•or* operator 3 



POLY array function 
PROD array function 
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SQRT 

complex arguments 35 
real arguments 19,20 

SUB SIR 

bit string 5 
character string 6 

SUM #5 



range of arguments 

in arithmetic operations and 
functions 8 

in array indexes and functions 

in mathematical functions 17 
relative error 19 
REPEi^T 

bit string 4 

character string 6 



i*3 



TMi 

complex arguments 4 

real arguments 25,26 
TAND (real arguments) 25,26 
TANH 

complex arguments 4 

real arguments 30,31 
truncation 

in array indexes and functions 43 

in math^natical functions 18 
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